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Abstract. The aim of this paper is to provide a handy tool to compute the impact of Type la SN (SNIa) events on the evolution 
of stellar systems. An effective formalism to couple the rate of SNIa explosions from a single burst of star formation and the 
star formation history is presented, which rests upon the definition of the realization probability of the SNIa event (Ai a ) and the 
distribution function of the delay times (/i a (r)). It is shown that the current SNIa rate in late type galaxies constrains Ai a to be 
on the order of 1CT 3 (i.e. 1 SNIa every 1000 M Q of gas turned into stars), while the comparison of the current rates in early and 

■ late type galaxies implies that fi d ought to be more populated at short delays. The paper presents analytical formulations for the 
description of the fj a function for the most popular models of SNIa progenitors, namely Single Degenerates (Chandrasekhar and 

. Sub-Chandrasekhar exploders), and Double Degenerates. These formulations follow entirely from general considerations on 

^ ' the evolutionary behavior of stars in binary systems, modulo a schematization of the outcome of the phases of mass exchange, 

\& | and compare well with the results of population synthesis codes, for the same choice of parameters. The derivation presented 

here offers an immediate astrophysical interpretation of the shape of the /i a functions, and have a built in parametrization of the 
key properties of the alternative candidates. The important parameters appear to be the minimum and maximum masses of the 
components of the binary systems giving rise to a SNIa explosions, the distribution of the primary mass and of the mass ratios 
in these systems, the distribution of the separations of the DD systems at their birth. The various models for the progenitors 
correspond to markedly different impact on the large scales; correspondingly, the model for the progenitor can be constrained 
by examining the relevant observations. Among these, the paper concentrates on the trend of the current SNIa rate with parent 
' galaxy type. The recent data by Mannucci et al. 120051 favor the DD channel over the SD one, which tends to predict a too steep 

distribution function of the delay times. The SD scenario can be reconciled with the observations only if the distribution of the 
mass ratios in the primordial binaries is flat and the accretion efficiency onto the WD is close to 100%. The various models 
are characterized by different timescales for the Fe release from a single burst stellar population. In particular the delay time 
| within which half of the SNIa events from such a population have occurred, ranges between 0.3 and 3 Gyr, for a wide variety 

■ of hypothesis on the progenitors. 
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1 . Introduction recorded in its stars. This argument is at the basis of the eval- 
uation of the formation timescales of stellar systems from the 

The evolution of the rate of Type Ia Supernovae (SNIa) with . „ , , . , , , , . .. . 

J * r * a to Fe abundance ratios, and has been used to estimate the 

time is a fundamental ingredient for the study of a variety of as- £ ir.uuir i r j--.ii 

° , . J formation timescale of the halo of our galaxy from individual 

trophysical issues, ranging from the chemical evolution of stel- . „ , , , . „ . . . „ „ ■ i lri0/; i „ 

r J ' & & stellar abundances (e.g. Matteucci & Greggio 1986 1, as well as 

lar systems, to the interpretation of the SNIa rates as a function . . c , , c .. , c c t ,, 

f d h'ft infer short formation timescales for Es from the Magnesium 

and Iron indices in their spectra (e.g. Matteucci 1994 Greggio 

Indeed, while type II Supernovae (SNII) precursors are [1997! Thomas et al.EOOll. SNIas are also thought to produce a 

short lived, massive stars, so that their rate evolves (almost) f .. . . • , . ■• • , 

' v ' major fraction of the iron in the intracluster medium in galaxy 

in pace with the rate of star formation, SNIa come from bi- dugterg ( Matteucci & Ve ttolani 1988; Renzini 1997J; the 

nary systems with a wide range of lifetimes, as indicated by temporal behayior of ^ SMa fate ^ ^ impact Qn ^ ^ 

their occurrence in both late and early type galaxies. Thus, abundance of the intracluster gas at high reds hift. 
the SNIa products are released to the interstellar medium over 

, .. , j . CXTTT . , c - ,u 1 Type Ia SNe provide the heating mechanism of the mass 

longer timescales, compared to SNII products. Since the a el- ,r F & 

j j u omtt u i cxtt • lost by stars in Ellipticals, an therefore determine the dynam- 

ements are mostly produced by SNII, while SNIa are impor- 1 r 1 

. . ... . t ■ , , c .. ■ , ical evolution of the gas in these galaxies. Following Renzini 

tant contributors of iron, the shorter the formation timescale . . b b b 



c , , • ,. , . — , , t . ( 1993), this evolution depends on the balance between the rate 

of a stellar system is, the higher the a to Fe abundance ratio ' ' r 

of mass return, and the secular evolution of the SNIa rate: if 

Send offprint requests to: L. Greggio the latter decreases faster than the former, the early stages are 
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characterized by supersonic winds, which then turn to subsonic 
outflows, and eventually to inflows (Ciotti et al. 1991). The re- 
verse sequence applies instead in the case of a mild secular evo- 
lution of the SNIa rate, with early inflows eventually turning 
into winds, such as in the models of Loewenstein & Mathews 

119871 . 

Finally, cosmological applications of the SNIa rate include 
the possibility of deriving clues on the SNIa progenitors, and/or 
constraining the cosmic star formation history from the evolu- 
tion of the SNIa rate with redshift (Madau et al. 1998 1 Strolger 
et al 2004). Understanding the evolutionary path which leads 
to a Type la explosion is of great importance to assess the use 
of SNIa as distance indicators, and to derive the cosmological 
parameters (Riess et al. 1998, Perlmutteret al. 1999). 

To address these issues and quantitatively interpret the re- 
lated observations we need a suitable description of the evo- 
lution of the SNIa rate from a single burst stellar population. 
From a theoretical point of view, this rate is difficult to derive, 
first because the nature of the progenitors of SNIa events is still 
an open question, and, second, because any theoretical rendi- 
tion is highly model dependent. While in the literature there is a 
general consensus that SNIa originate form the thermonuclear 
explosion of Carbon and Oxygen (CO) White Dwarfs (WD), 
various evolutionary paths may lead to such event. Common to 
all the models is the first part of the evolution, dealing with a 
close binary system with a primary component less massive 
than ~ 8 M Q , so that it evolves into a CO WD. When the 
secondary component expands and fills its Roche Lobe, the 
primary may or may not accrete the matter shed by its com- 
panion. If the accretion rate is approximately ~ 10~ 7 M Q /yr 
(Nomoto 1982), the accreted matter burns on top of the WD, 
the object remains confined within its Roche Lobe and grows in 
mass (Whelan & Iben[T373}. However, Iben & Tutukov {T984) 
pointed out that in most cases the secondary expands at such 
a high rate that the accretion rate exceeds the mentioned limit, 
implying that a common envelope (CE) forms around the two 
stars. Orbiting inside the CE the two cores spiral in, and orbital 
energy is transferred to the envelope which is eventually lost. 
The system emerges form the CE phase as a close double WD, 
which will merge due to the emission of gravitational wave ra- 
diation. Another interesting possibility is that, when subject to 
a large accretion rate, the WD develops a strong radiative wind, 
to the effect of stabilizing the mass transfer, thus allowing the 
WD to grow in mass, and eventually explode (Hachisu, Kato 
& Nomoto 1996). Within all scenarios explosion occurs either 
when the CO WD reaches the Chandrasekhar mass and Carbon 
deflagrates at the center (Chandra exploders), or when a mas- 
sive enough helium layer is accumulated on top of the CO WD, 
so that Helium detonates, inducing off center Carbon detona- 
tion (e.g. Woosley & Weaver [l994i before the Chandrasekhar 
mass is reached (Sub-Chandra exploders). 

Different arguments can be found in favor or against both 
scenarios (e.g. Livio|2p01), generally referred to as Single 
Degenerate (SD) and Double Degenerate (DD), depending 
on whether the SNIa precursor is a system with one or two 
WDs. Briefly, the SD model is supported by the observa- 
tional detection of several classes of objects that can be con- 
sidered as potential SNIa precursors of the SD variety, i.e. 



Cataclysmic Variables, Symbiotic Stars, and Supersoft X-Ray 
Sources (Munari & Renzini [79921 Kenyon et al. [T9931 Van 
den Heuvel et al. 119921 Rappaport, Di Stefano & Smith 
1994). Additional support to the Cataclysmic Binaries chan- 
nel came recently from the detection of a candidate companion 
to Tycho's supernova (Ruiz-Lapuente et al. 2004). On the other 
hand, the fine tuning of the mass accretion rate limits consider- 
ably the volume in the parameter space for a successful SNIa 
explosion in the SD model. As a consequence, it seems likely 
that only a small fraction of events can be realized through this 
channels in our galaxy (Fedorova, Tutukov & Yungelson 2004; 
Han & Podsiadlowski 2004, but see Hachisu, Kato & Nomoto 
1999 for a different point of view). 

Several attempts have been made to establish the binary fre- 
quency among White Dwarfs, and to determine the distribu- 
tion of total masses and periods of the binary systems, in order 
to assess the likelihood of the DD channel as SNIa progenitor 
(Robinson & Shafter [T987l Bragaglia et al. [T990l Foss, Wade & 
Green fT99ll Saffer, Livio & Yungelson fT998l Maxted & Marsh 
1999 ). To date, the most comprehensive effort to find SNIa pre- 
cursors among DD systems is the SPY project (Napiwotzki 
et al. 2001), whose results have been recently summarized in 
Napiwotzki et al. (2004): many close DD systems have been 
found, with one very good candidate SNIa precursor, with a 
total mass exceeding the Chandrasekhar limit and expected to 
merge within a Hubble time. In addition, a few other systems 
come close to qualify as SNIa precursors. In general, it seems 
that the masses and period distributions of the binary WDs 
confirm the prediction of the population synthesis models; ac- 
cording to these models the DD evolutionary channel is able to 
provide enough merging events to match the current SNIa rate 
measured in our galaxy, which is similar to the typical SNIa 
rate in Spirals. On the other hand, theoretical calculations show 
that the merging of two massive WDs may lead to accretion in- 
duced collapse, rather than to SNIa explosion (Saio & Nomoto 
1998), so that the ultimate fate of these candidates may be a 
neutron star. 

The various models for the progenitors, SD or DD, un- 
dergoing Chandra or Sub-Chandra explosions, correspond to 
rather different temporal behavior of the SNIa rate (see e.g. 
Fig. 2 in Yungelson & Livio 2000 hereafter YL). In the cur- 
rent literature there are several examples of theoretical com- 
putations of the SNIa rate performed with population syn- 
thesis codes: starting from a primordial distributions of bi- 
nary masses, mass ratios and separations, the computations fol- 
low the evolution of the stellar systems under some prescrip- 
tions for the mass exchange between the binary components, 
to determine the final outcome (Tutukov & Yungelson 1994; 
Yungelson et al.Q594; Ruiz-Lapuente, Burkert & Canal 1995, 
Han, Podsiadlowski & Eggleton [T995J Nelemans et al. l200Tl 
De Donder & Vanbeveren 2003 ). The results of these simu- 
lations depend on a variety of input parameters and assump- 
tions, whose role is difficult to gauge, so that they are not suited 
to easily explore the parameter space for the SNIa progeni- 
tors' models. In addition, the population synthesis codes yield 
numerical outputs, which are difficult to incorporate in codes 
which follow the evolution of galaxies or of galaxy clusters. 
Indeed, most of the astrophysical applications of the SNIa rates 
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in the literature are based either on the analytical formulation 
by Greggio & Renzini ( 1983}, or on the parametrization pro- 
posed by Madau et al. ( 1998 1. However, the former is derived 
only in the framework of the SD model; the latter is a conve- 
nient mathematical expression, but it is only marginally related 
to the physics of stellar evolution. 

In this paper I provide relatively simple analytical formula- 
tion for the SNIa rate, which allows us to identify the most crit- 
ical parameters and should help restrict the choice among the 
candidate precursors. Both the SD (Chandra and Sub-Chandra) 
and the DD (only Chandra) models are considered, so as to pro- 
vide a handy way for investigating on the impact of the different 
SNIa models on astrophysical issues for which the SNIa rate is 
important. Section 2 presents a coherent formalism to couple 
a particular SNIa model with the star formation history of a 
system. The analytic expressions for the SNIa rate for the SD 
and DD models are derived in Sections 3 and 4 respectively. 
Readers mostly interested in the main results may skip these 
sections, since a general description of the analytic /i a function 
appears at the beginning of Section 5, where they are compared 
to the predictions of population synthesis codes. In addition, 
Section 5 presents an attempt to constrain the SNIa progeni- 
tors from the systematic trend of the SNIa rates with galaxy 
type. Finally, some concluding remarks appear in Section 6. 
The mathematics used to derive of the analytic relation for the 
DD model is (mostly) described in the Appendix, for an easier 
readability of the text. 

2. Formalism 

A convenient formulation of the SNIa rate to follow the evolu- 
tion of stellar systems rests upon the definition of the distribu- 
tion function of the delay times, i.e. the time elapsed between 
the birth of a SNIa progenitor and its explosion. I indicate this 
function with /i a (r), defined in the range (t;, t x ) , respectively 
the minimum and maximum possible delay times, and consider 
/ia(i") normalized to 1 : J t 1 fjj) ■ dr — 1 . The minimum de- 
lay time Tj is the minimum evolutionary lifetime of the SNIa 
precursors: for the SD model this is the nuclear lifetime of the 
most massive stars which produce a WD, that is an ~ 8 M star, 
which evolves off the MS in ~ 0.04 Gyr. For the DD model, 
Tj could be appreciably larger than this because of the addi- 
tional gravitational delay, i.e. the time taken by the DD system 
to merge due to the gravitational wave radiation. The maxi- 
mum delay time r x is quite sensitive to the model for the SNIa 
precursor, as will be seen later. At this point I just notice that, 
if elliptical galaxies formed the bulk of their stars in a short 
initial burst, the maximum delay time of their inhabiting SNIa 
precursors must be on the order of a Hubble time, or more. 

At a given epoch f, the contribution to the SNIa rate from 
progenitors with delay times in the range (r, r+dr) is 

dh Ia = hjt-r)x A Ia (f - r) x fijj) dr (1) 

where n*(f - r) is the birth rate at epoch (f - r), and Ai a is 
the realization probability of the SNIa scenario from the stellar 
generation born at (f - t) 1 . In Eq. {1} I have considered the 

1 If N t is the number of stars born at epoch (t-r), Ai a x N* is the 
number of SNIa events produced by this stellar generation ever. 



possibility of variations of Ar a during the galaxy evolution. For 
example, one could expect larger realization probabilities of the 
SNIa channel at higher metallicities because, during their evo- 
lution, stars expand to larger radii (Greggio & Renzini 1990 1, 
hence have a better chance to fill their Roche Lobe. In addi- 
tion, at high metallicity the accretion onto the pre-supernova 
WD could be more efficient (Hachisu et al. U996> . also imply- 
ing larger values for Ai a . 

Following Tinsley ( 1980 1 notation: 

h+(t -T) = iff(t-r)X I <f>(m) dm (2) 

where <p{m) is the initial mass function (IMF) by number, and 
;«, and m s are the lower and upper mass limits. As usual, I adopt 
a power law IMF, with total mass normalized to 1, so that the 
star formation rate (SFR) iff is the mass that goes into stars per 
unit time. It follows: 

H*(? - T) = tff{t -T)Xk a (3) 

where k a , which is the number of stars per unit mass in one 
stellar generation, depends on the IMF. For example, for a mass 
distribution 0(m) oc m~ a ranging from 0.1 to 120 M , k a is 
equal to 2.83 and 1.55, respectively for Salpeter (a = 2.35) 
and Kroupa J2"0"0"T1 > (a = 2.3 in m > O.5M and a = 1.3 in 
m < O.5M ) IMFs. 

By substituting Eq. Q into Eq. Q, and summing over all 
the contributions from the past stellar generations, the SNIa 
rate at epoch t is: 

4>(t-r) A Ia (f-r) / Ia (r) dr. (4) 

This equation allows one to insert consistently the SNIa 
events in codes which describe the evolution of galaxies, once 
the distribution of the delay times fi a and the fraction Ai a are 
specified. Notice that these are the only results of the model- 
ing of binary populations of SNIa precursors which impact on 
the evolution of stellar systems. In other words, for the astro- 
physical applications, the particular prescriptions used in the 
population synthesis codes are of limited interest, while most 
important is the distribution function of the delay times and 
the total realization probability of the SNIa scenario out of one 
stellar generation. 

Eq. © can be easily specified for a single generation of 
stars by considering a star formation episode started at t — 
and proceeded at a constant rate ipg for a very short time (Afe). 
In this case, the integrand function is non zero only in a narrow 
age range around r = t so that: 

n,Jt) = k a -ifr B - Af B • A Ia (f = 0) ■ / Ia (r = f) 

= k a ■ At B ■ Aia • Mi) (5) 

where Mb is the mass that went into stars in the burst. Thus 
the SNIa rate following an instantaneous burst of star forma- 
tion is proportional to the distribution function of the delay 
times through a factor which is the product of the total stel- 
lar mass formed in the burst, times the realization probability 
of the SNIa scenario, times the number of stars per unit mass, 
characteristic of the IMF. 
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2.1. Observational constraints 

Eq.0 shows that the SNIa rate results from the convolution of 
the distribution function of the delay times and the SF history: 
in order to derive clues on the former from the observed rates 
I consider a family of models, characterized by a constant SFR 
(t^o) starting at t = 0, and with variable durations (At). I also 
assume that early (late) type galaxies are represented by models 
with short (long) Af. Although the SF history in real galaxies is 
much more complex than that, this schematic description leads 
to interesting relations between the observed SNIa rates and 
the key quantities Ar a and /r a , which will not be far from those 
that one can obtain with a thorough modeling of the galaxy 
evolution. 

For this family of models, the minimum delay time which 
contributes to the rate at epoch t is either T\ (i.e. the absolute 
minimum delay of the SNIa progenitors), or the time elapsed 
from the end of the burst: t - At. Taking this into account, and 
neglecting the temporal dependence of the factor Ai a , Eq. © 
becomes: 



hia(t) =ka-Au-ifro- I /ia(r) dr. 

Jmax(Ti,/-Ar) 



(6) 



The two options for the lower integration limit correspond 
to two different regimes: 

1. when t; > t- At , i.e. t < t; + Af, we basically map the SNIa 
rate while the SF is still ongoing; 

2. when t\ < t - At, i.e. t > t; + Af, we get the SNIa rate after 
the burst is completed. 

The transition between regime 1) and 2) occurs at older ages 
when the SF episode lasts longer (i.e. when Af is larger). For a 
currently ongoing SF the transition has not yet occurred. 

In regime 1), relevant for late type galaxies, the SNIa rate 
is given by: 



nf a (t) 



k a ■ A Ia • if/ 



■min(/',T x ) 



ZiaWdr 



(7) 



showing that it increases with time, as systems with in- 
creasingly longer delays contribute to the explosions, up to 
f = r x . From then on the SNIa rate stays constant, and equal 
to 



h L la (t) = k a ■ Ai a • iAo 



(8) 



since /i a is normalized to 1 . Notice that this expression is valid 
only when Af > r x , that is when the SF episode lasts long 
enough to include all the possible delay times. In this case the 
current SNIa rate gives indications on the realization probabil- 
ity A Ia . For example, Cappellaro, Evans & Turatto ( 1999 ) report 
an observed SNIa rate in late type galaxies of 0.2 SNu 2 , i.e. 
0.2 • (Lb/^b,©) ■ 10~ 12 events per year. Then, approximating fa 
with the ratio between the galaxy's stellar mass (M) and its age 
(f) Eq. © yields: 



fft _ 3 0.2 (M 



X f G y 



(9) 



with Lb in solar units. Notice that the mass M in Eq. corre- 
sponds to the integrated star formation rate, and then it differs 
from the actual galaxy stellar mass because of the gas returned 
by the stars as they evolve. Using Maraston ( 1998 2005) sim- 
ple stellar population models with solar metallicity, I find that 
for a constant SFR the Ai/Ls ratio at an age of 12 Gyr is ~ 2 or 
1.5, respectively for a Salpeter and Kroupa IMF. By inserting 
these values into Eq. l|9} coupled with the appropriate k a , I get 
Ai a ~ 5 x 10~ 4 or 10~ 3 , respectively for Salpeter and Kroupa 
IMF. Thus, one SNIa event is produced for every 2000 (1000) 
M Q of stars formed with a Salpeter (Kroupa) initial mass dis- 
tribution. Although this is only a rough estimate, it appears that 
the current rates measured in late type systems require rela- 
tively small realization probabilities of the SNIa channel. For 
comparison, the realization probability of SNII, estimated as 
the fraction of stars more massive than 8 M in a stellar popu- 
lation, is ~ 0.003, 0.007 respectively for Salpeter and Kroupa 
IMF. Therefore, out of a stellar population of 1000 M Q roughly 
1 SNIa, and 7 SNII are produced. 

Considering now regime 2), relevant for early type galaxies, 
as long as f < r x Eq. (|6jl yields: 



hf a (t) = k a ■ A Ia ■ iff ■ Af • (/ia> f -A»,( 
= k a ■ A u ■ M ■ (fu)t-At,t 



(10) 



where (/i a )r-A/,/ is the average of the distribution function of the 
delay times in the range (f - Af; f). 

After the burst is completed, the SNIa rate scales according 
to the /ia function, and goes to zero at f = r x + Af. Early type 
galaxies are currently in this age regime, so that their SNIa rate 
yields information on the value of the distribution function of 
the delay times at ages in the vicinity of the galactic age. This 
constraint is better illustrated by constructing the ratio of the 
SNIa rate in early and late type systems at the current epoch. 
Dividing Eq. fllOi by Eq. 0, and considering the rates in SNu, 
we get: 



„ _ h Ia.SNu (M/L B ) E </la),-A,,, 

SNu = " (M/i*r x </iak, ' 



(id 



having assumed that the IMF and the Ai a factors are the same 
in early and late type galaxies. According to Cappellaro et al. 
( 1999), the current epoch SNIa rate in early and late type sys- 
tems is the same. This puts a constraint on the ratio between 
the average values of the function fi d : 



{fu)t-h 



(M/L B ) L 



2 1 SNu is one event per century per 10 10 L BG . 



</la>T„, (AI/Lb) 

Maraston ( 1998 2005 ) models of simple stellar populations 
at an age of 12 Gyr and solar metallicity have Ai/Lg~ 13,10 
respectively for Salpeter and Kroupa IMF, again having con- 
sidered the mass initially transformed into stars 3 . Combining 

3 Notice that at 12 Gyr, ~ 30 % of the initial mass has been returned 
to the interstellar medium by the evolving stars, for a Salpeter IMF. 
The return fraction is instead 45 % if Kroupa IMF applies. Therefore, 
the ratio between the current stellar mass and the blue luminosity pre- 
dicted by the models is ~ 9, 5 for the two IMFs. 
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this estimate of the M/Lb ratio for early type systems with the 
values quoted above for the late type systems, it turns out that 
the current SNIa rates indicate: 

(/la)(-A/,/ „ _ 

Thus the value of /j a at late delay times is substantially 
smaller than its average value over the whole range f); in 
other words, the distribution function of the delay times de- 
creases with time. This means that the majority of SNIa pre- 
cursors are relatively short lived. 

It is worth emphasizing that if the distribution of the delay 
times were flat (i.e. if young and old systems were equally effi- 
cient in producing SNIa events), the observed ratio ^snu would 
be of the order of 5-10, given the ratio of the Ai/Ls values 
for early and late type galaxies. Therefore, the fact that the /i a 
function must be decreasing with increasing delay time is a 
very robust conclusion. 

2.2. Stellar evolution predictions 

The distribution function of the delay times and the realization 
probability of the SNIa channel from a stellar generation can be 
derived from the theory of the evolution of binary systems. As 
mentioned in the Introduction, in the literature there are sev- 
eral examples of the SNIa rate in stellar systems computed 
with population synthesis techniques (see Yungelson 2004 1. 
Typically, these models predict a /i a function characterized by 
an early maximum, and a late epoch decline, while the realiza- 
tion probability is indeed on the order of 10~ 3 . However, these 
results depend on the adopted parameters of the simulations, 
like the star formation history, the distribution functions of the 
separations and masses of the primordial binaries, and on the 
specific prescriptions for the evolution during the hydrodynam- 
ical phases of the mass transfer. Therefore, (i) the resulting /i a 
functions are model dependent; (ii) the role of the various input 
parameters on the output is far from straightforward. 

On the other hand, on general grounds, Ciotti et al. (I199U 
give a motivation for the late epoch decline related to the tem- 
poral behavior of the clock of the explosions. Indeed, some 
general considerations can be made which strongly character- 
ize the shape of the distribution function of the delay times, as 
I show in the following sections. I consider separately the two 
main categories of SD and DD progenitors, and derive analyti- 
cal formulations for the /i a (r) function, in the attempt of clari- 
fying the role of the important stellar evolution parameters. 

3. Single Degenerates 

In the SD model, the clock of the supernova event is set by 
the evolutionary lifetime of the secondary. A fit to Girardi et 
al. (|2000 1 solar metallicity tracks yields the following relation 
between stellar mass (in M ) and Main Sequence (MS) lifetime 
(tms, in years ), valid in the range (0.8 < m 2 < 8)M Q , which 
corresponds to (0.04 < r M s 5 25) Gyr : 

log m 2 = 0.0471 x (log r MS ) 2 - 1-2 x log r MS + 7.3. (12) 




log T 

Fig. 1. Evolutionary mass (solid line) and its derivative (dot- 
dashed line) as a function of the MS lifetime (in years), for 
solar metallicity tracks. The derivative (to be read on the right 
axis) drops by ~ 4 orders of magnitude as the MS lifetime in- 
creases from ~ 40 Myr to ~ 20 Gyr. The dashed line represents 
the distribution function of the delay times if the secondaries 
follow a Sapeter IMF (see text). 

The total delay time is basically equal to the MS lifetime of 
the secondary component of the binary system, the Post-MS 
phase being in any case much shorter than the hydrogen burn- 
ing lifetime. Therefore, for one stellar generation, the number 
of explosions within (r, r+dr) is proportional to the number 
of binaries with secondary mass between m 2 and m 2 +dm 2 such 
that the evolutionary lifetime of m 2 is t: 

/if (t) |dr| oc n(m 2 ) \dm 2 \, (13) 

which implies: 

/ t f (t) oc n(m 2 ) ■ \m 2 \. (14) 

Figure n snows the evolutionary mass from Eq. i ll 2t and its 
derivative as a function of the delay time t- tms- The deriva- 
tive is very well approximated by the power law \m 2 \ oc T - iM 
over the whole range from 40 Myr to > 15 Gyr. In Eq. (I14> . 
as t increases, the factor \m 2 \ strongly decreases. However, 
longer delay times correspond to smaller evolutionary masses 
m 2 , and the final shape of the distribution function depends on 
the distribution of secondary masses as well. The dashed line 
in Fig^shows the result obtained with a Salpeter distribution 
n(m 2 ) oc OTj 2 35 , clearly not steep enough to imply a / Ia D (T) in- 
creasing at late delay times 4 . 

However, while the distribution function of the secondary 
masses of primordial binaries could well be represented by a 

4 Since 1112 oc t~ 03 at ages older than 1 Gyr, mass distributions as 
steep as n(m 2 ) oc m^ 4 ' 5 are required in order to have a flat /^ D at late t. 
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Salpeter IMF, the distribution n(mi) in Eq. (I14> refers to the 
secondaries in systems which eventually give rise to a SNIa 
event, and as such suffers from some limitations. This is illus- 
trated in the next section. 

3.1. The distribution function of the secondary masses 

The commonly adopted scheme in the population synthesis 
computations deals with binary systems in which the primaries 
follow a power law distribution with slope -a, and the mass 
ratios q = m 2 jm\ follow a power law distribution with a 
slope y. Then, the number of binaries with primary masses in 
(mi; mi + dm\) and secondaries in (m 2 \ m 2 + dm 2 ) is : 

n(mi, m2) dmi d;«2 = n(m\) f{q)Am\Aq. (15) 

The distribution function of the secondaries in SNIa pro- 
genitor systems, is obtained by summing over all possible pri- 
maries, ranging from a minimum value (mi ;) to 8 M : 

f 8 - a+ +1 

n(m 2 ) °c I m^flJj r d»ii 

cc m 2 ° ■ [(m 2 /m ui r + ? - (m 2 /&T + ?]. (16) 

The minimum mass for the primary (m^) is either equal to 
1112 (as in any binary system), or more massive than this, if the 
primary has to produce a sufficiently massive CO WD in order 
to lead to a SNIa explosion. Thus: 

m\i — max(m2,mi in ) 

where m\ R is the mass of the primary following a specific con- 
straint on the minimum mass of the CO WD. The restriction 
of the integration to systems with primaries more massive than 
mi n has important consequences on the distribution function of 
the delay times, and then, on the SNIa rate past a burst of SF, 
as will be seen in the following. For Chandrasekhar explosions, 
mi itl is derived by requiring that 

™WD,n + £ ■ '«2,e = 1.4 (17) 

where ;«wD,n is the minimum acceptable mass for the WD, m2, e 
is the envelope mass of the evolving secondary, and e is an effi- 
ciency parameter. In principle mz, e varies with the evolutionary 
stage at which the second Roche Lobe Overflow (RLO) oc- 
curs, and then depends on the separation of the binary system. 
In practice, since the second mass transfer will not necessar- 
ily occur with 100% efficiency, this detail can be neglected, to 
consider a representative relation between m2 and m 2& . For ex- 
ample, the results of case B RLO in Nelemans et al. (2001) 
can be represented by the following relation between the initial 
mass (m2) and its remnant (m 2fi ): 

m 2 ,c = max{0.3; 0.3 + 0.1(m 2 - 2); 0.5 + 0.15(m 2 - 4)}. (18) 

By considering m 2t = m 2 - '«2,c it turns out that smaller secon- 
daries typically have smaller envelopes to donate to their WD 
companion. As a consequence, the minimum acceptable WD 
mass m WDin increases as m 2 decreases, and so does the mass of 



its progenitor mi A . An additional lower limit on the mass of the 
primary comes from the requirement that its remnant should 
be a CO, rather then a Helium, WD. Stars less massive than 
~ 2 M Q develop a degenerate Helium core, and evolve along 
the Red Giant Branch up to the Helium flash, which occurs 
when the stellar radius is of a few hundred R Q (see e.g. Fig. 1 in 
Yungelson 2004). Therefore, these stars can provide CO WDs 
only if the separation of the primordial binary exceeds ~ 400 
R . Since the distribution function of the primordial separa- 
tions Ao scales as Ag 1 (see, e.g. Iben & Tutukov 1984) systems 
with a primary less massive than 2 M are much more likely 
to produce a Helium rather than a CO WD. For this reason, the 
contribution to SNIa from systems with m\ smaller than 2 M 
is neglected here, and I consider: 

mi >n = max{2.,2. + 10.(m W D,n - 0.6)} (19) 

where the relation between m\ fi and its remnant mwD.n de- 
scribes the results of case C RLO in Nelemans et al.(2001 1, and 
represents very well the empirical determination of the initial - 
final mass relation by Williams, Bolte & Koester (2004). 

In summary, as m2 decreases, the lower limit of integra- 
tion in Ea. (ll6> is first set to m\j=m 2 , down to m2=2 M Q ; then 
it stays constant and equal to 2 M Q until the minimum WD 
mass for a successful SNIa event (mwD.n from Ea.( ll7» be- 
comes larger than 0.6 M . From that point on, m\\ increases 
with decreasing m 2 . At some value of m 2 , m\j becomes equal 
to 8M : this marks the minimum secondary mass suitable for a 
successful SNIa event. 

The scheme adopted in Greggio and Renzini (1983} is 
slightly different from the one just described. It assumes that 
the total mass of the primordial binary Mb follows a power law 
with slope -a, and that the mass ratio fi-m 2 /M\, is distributed 
according to f(p) cc yu y . In this case: 

n(m 2 ) cc m r 2 M b (a+y+1) dM b 

JM bi 

cc m 2 a ■ [(2m 2 /M^) a+y - (2m 2 /M b ^] (20) 

where Mb, s = 8+m2 and M b j=mn+m 2 . 

It is worth remarking that in the original formulation by 
Greggio & Renzini d!983> . Mb,; was required to be larger than 
a minimum value (e.g. 3 M ) irrespectively of the mass of the 
evolving secondary. This fixed limit does not describe the one 
to one correspondence between the mass of the evolving sec- 
ondary and delay time of the SNIa precursor. At given delay 
time, the amount of mass that can be donated to the accreting 
white dwarf is virtually fixed, which implies a minimum white 
dwarf mass, in order to add up to the Chandrasekhar limit. A 
revision of the Greggio & Renzini ( 1983) SNIa rate has been 
presented in Greggio ( 1996 ). 

Equations (1161 and J20i are very similar, but not quite the 
same. In both cases, the distribution function of the secondaries 
follows a power law with slope -a, corrected by a factor which 
describes the limitations imposed on the masses of the pri- 
maries in order to secure the SNIa explosion. As t increases, 
m 2 decreases and, as long as m\;{= m 2 , the correction factor 
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Fig. 2. The distribution function of the secondary masses in the 
SD model, for Chandra exploders as obtained with Eqs. dl6> 
and (I20> for selected values of the a and y parameters. Solid 
lines have been obtained with e = 1, dashed lines with e = 0.5 
in Eq. dl7> . The black curve shows the derivative of on a 
linear scale (right axis). 

increases toward unity. As soon as my becomes greater than 
m.2, the correction factor starts decreasing, due to the loss of 
systems with primary mass between m.2 and my. 

Figure (0 shows the distribution (in arbitrary units) of the 
secondaries in systems which lead to Chandrasekhar explo- 
sions, under a variety of hypothesis for the power law slopes 
a and y, and for the two formulations given by Eqs. (II 61 and 
J20i . as labeled. 

The distribution function of the secondaries in systems 
which produce a SNIa is shaped after the behavior of the min- 
imum mass for the primary as mi decreases: it shows a first 
abrupt change of the slope when m-i drops below 2 M , that is 
when my is set to 2 M independently of mo. A second abrupt 
change of the slope appears when the envelope mass of the sec- 
ondary is so small that my > 2M is required in order to build 
up to the Chandrasekhar mass. From this point on, the distri- 
bution function of the secondaries steeply decreases as my in- 
creases. Smaller values for e imply an earlier occurrence of this 
regime, which in the illustrated case (e=0.5) appears soon after 
mi has gone below 2 M . 

For the same values of a and y, the Greggio & Renzini 
( 1983 1 formulation (red curves) yields a larger fraction of sys- 
tems with low secondary mass, compared to the formulation 
generally adopted in the population synthesis codes. As a re- 
sult it leads to comparatively larger rates at late epochs. 

The distribution n(mi) also appears very sensitive to the a 
and y parameters. Most noticeably, the flatter y, the larger the 
fraction of systems with small secondaries. Figure[2]also shows 
the time derivative of mi: the combination of the two factors 



Fig. 3. Distribution function of the delay times for the SD 
model for the same choices of the parameters as in Fig. |2] In 
addition, the dot-dashed lines show the distribution function for 
the Sub-Chandra models in the case (a, y)=(2.35,l), using Eq. 
(O (blue) and Eq. (|20j (red). 



clearly produces an early maximum for the /? function which 
is given by Eq.( ll4t . 

3.2. The distribution function of the delay times 

Figure [3] shows the resulting distribution function of the delay 
times, with the same color and line coding as in Fig. 13 The two 
abrupt changes of the slope just reflect those appearing in the 
n(m2) function. Figure [3] also shows the distribution functions 
of the delay times expected for the Sub-Chandra exploders as 
dot-dashed lines. These are obtained from Eqs d!6t and i2Q\ 
with the following criterion for my: a minimum WD mass of 
0.7 M Q is required to secure the explosion, along the lines sug- 
gested by Woosley & Weaver ( 1994 1 models. This corresponds 
to a minimum primary mass of 3 M (see the second term of 
the RHS in Eq. dl9t ). so that the first cusp occurs as early as 0.4 
Gyr. In order to have a Sub-Chandra SNIa, a Helium layer of 
about 0. 15 M needs to be accumulated on top of the CO WD, 
which corresponds to (e • m2, e ) > O.15M . This is a limit on 
mi which implies that the single burst SNIa rate for this model 
drops to zero at a delay time equal to the evolutionary lifetime 
of such n%2. Since the evolutionary mass at 15 Gyr is ~ 0.9 M , 
this limit is inactive for the whole Hubble time, provided that 
e ^ 0.3, as is considered for Fig.|3 

The distribution function of the delay times for the SD 
model is basically shaped according to the limits on the mass 
of the primaries in systems which eventually produce a SNIa. 
Four regimes can be identified: 
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1) Regime A, in which shows an initial steep rise followed 
by a mild decrease. In this regime m\$=m,2, implying that at 
longer delay times there is an increasing range of ni\ which 
contribute to the SNIa explosions. The (wide) maximum results 
from the interplay between n(m 2 ) and \m-i\, the latter being a 
decreasing function of t. 

2) Regime B, which sets in when nij falls below the minimum 
primary mass suitable for the explosion. Such minimum mass 
is chosen here equal to 2 and 3 M Q respectively for Chandra and 
Sub-Chandra exploders. In Regime B, mij is kept constant and 
equal to this minimum primary mass; the distribution function 
of the delay times decreases with a steeper slope compared to 
Regime A, due to the loss of systems with ra 2 < raj < 2. It 
can be shown that in this regime, the function is mostly 
sensitive to the parameter y, according to log ff® oc T -i-44-o.3y 

3) Regime C, starting when primaries more massive than 2 M Q 
are required in order to build up to the Chandrasekhar mass. In 
Regime C, mij increases from 2 to 8 M as r increases; conse- 
quently f^f shows a steep decline. For Sub-Chandra exploders 
Regime C does not exist. 

4) Regime D, in which /^ D =0: this sets in when mi t i=& M e is 
required to produce a Chandrasekhar explosion, or when the 
envelope mass of m 2 drops below the (0.15/e) limit. Regime D 
is not shown in Fig. [3] The value of r at which this regime sets 
in defines r x for the SD model. 

For the Chandra exploders, the delay time at which Regime 
C sets in varies with the efficiency of the accretion process. The 
response of the CO WD to accretion crucially depends on the 
accretion rate, as illustrated in Hachisu & Kato (2001 1. Only if 
the accretion rate is in a small range around 10 -7 M Q /yr does 
steady burning occur on top of the WD, so that the mass of 
the donor is efficiently used to increase the mass of the WD. 
For both lower and higher accretion rates, part of the donated 
mass is lost, either following H-shell flashes (novae explosion), 
for low accretion rates, or in a thick wind, for high accretion 
rates. Therefore, the accretion efficiency, i.e. the ratio between 
the accreted and the donated matter, is unlikely to be close to 
unity, as adopted for the solid lines in Fig. [3] The lower the 
accretion efficiency, the earlier regime C sets in, and regime B 
is suppressed when e < 0.5, as shown by the dashed lines in 
Fig El 

The shape of the /j^ D function derived here is quite general, 
since it just reflects the product of the time derivative of the 
evolutionary secondary mass, and the distribution w(m 2 ). 



4. Double Degenerates 

In the DD model, the first part of the close binary evolution 
is the same as in the SD model, but, following the expansion 
of the secondary component, a Common Envelope phase (CE) 
sets in, eventually leading to the complete ejection of the CE 
itself. In this scenario, accretion on top of the WD, if any, is 
neglected. The system emerges from the CE as a close binary 
composed of two WDs, which are bound to eventually merge 
due to the emission of gravitational wave radiation. If the to- 



tal mass of the DD system exceeds the Chandrasekhar mass, 
explosion occurs as a SNIa event. The delay time in this case is 



where r n is the nuclear lifetime of the secondary, and r gw is the 
gravitational delay (Landau & Lifshitz 1962 1: 



0.15 A 4 



(mm + OT 2 R)OTi R m 2 R 



Gyr 



(21) 



A, miR and m 2 R being respectively the separation and compo- 
nent masses of the DD system, in solar units. 

To derive the distribution function of the delay times one 
needs to map the distribution of the primordial systems in the 
space (Ao, rai,ra 2 ) into that of the final systems in the space 
(A, wiir, wi 2 r). Rather then performing Montecarlo simulations 
I consider here some general aspects, with the aim of charac- 
terizing the shape of the fi d function. 

I restrict to systems with 2 < m.\,m,2 < 8, from which 
most double CO WDs form, as argued in Sect. 3.1. Typically, 
the WD mass of both components ranges between 0.6 and 1.2 
M Q (see e.g. Eq. dl9». so that r n ranges between 0.04 and 1 
Gyr. The gravitational delay spans a large range, depending on 
the final separation: e.g., for a (0.7+0.7) M DD, r gw increases 
from 0.014 to 18 Gyr when A goes from 0.5 to 3 R e . The dis- 
tribution function of the delay times will depend on the distri- 
butions of both T n and r gw , at least up to total delays of ~ few 
Gyr, with early explosions provided by systems with short T n 
and r gw . Since the maximum nuclear delay does not exceed ~ 
1 Gyr, at long r the SNIa events will come from systems with 
long gravitational delays, that is DDs with wide separations and 
low masses. 



4.1. The interplay ofr a andr, 



gw 



The following useful approximation for the gravitational delay 
is justified in Appendix IA.1I 



= 0.6- 



Mi 



Gyr 



(22) 



where Mdd is the total mass of the DD system (again in solar 
units). This formula results from considering the restrictions on 
the masses of the two WD components. In this approximation 
the gravitational delay does not explicitly depend on ;m 2 r. Since 
T n only depends on m 2 , which is in tight correspondence with 
to 2 r, the relation between the two timescales is very weak: ba- 
sically, primordial systems with secondary mass ra 2 will evolve 
into a family of SNIa precursors spanning a wide range in sep- 
arations and total binary mass, that imply a wide range of grav- 
itational delays. 

The integrated distribution of the total delay times can now 
be constructed: the contribution (dn n ) of binaries with given r n 
to the number of systems with total delay time shorter than r 
is proportional to the fraction of them which have r gw shorter 
than (t - T n ). Indicating this quantity with g(j, r n ): 



dn n = «(t„) • g(r, t„) dr„ 



(23) 
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where n(r n ) is the distribution function of the nuclear 
timescales, proportional to the function described in Sect. 
3. The total number of systems with total delay time shorter 
than r is obtained by integrating the d« n contributions over the 
relevant r„ range: 



r p 

J n(r)dT = J 



min(r n x ,r) 



"(Tn) g(T, T n ) dr n 



where r n ,j and r n x bracket the range of nuclear timescales of 
the systems which end up in a successful SNIa event. 

Let's now indicate with n(r gw ) the distribution function of 
the gravitational delays of systems with a nuclear delay t„ (i.e. 
progeny of systems with a secondary mass ntj whose MS life- 
time is equal to r n ), and let r gWi ; and r gw x be the minimum and 
maximum gravitational delays of such systems. Formally, the 
fraction of systems which, having a nuclear delay equal to t„, 
have also a total delay shorter than r varies with r according to: 



g(j, T n ) 







■min(r gw . x ,r-r n ) 



"(Tgw) dr g 



for T< Tl 

for Tl < T< T2 

for T> T2 



(25) 



with Ti = T n + T gWi ; and T2 = T n + T gWjX . This equation merely 
expresses that (i) systems with nuclear delay r n have a total de- 
lay ranging between T n +T gWi ; and T n +T gWiX ; (ii) the total delay 
of such systems scales according to the distribution function of 
their gravitational delays. 

The shape of «(T gw ) will reflect the distribution of the fi- 
nal separations and of the DD masses. The interesting range 
of gravitational delays, which goes from ~ 0.01 Gyr (i.e. on 
the order of T n> ;) to over the Hubble time, is populated by low 
mass systems (M DD =1.4M ) with 0.5 < A/R Q < 2.8, and by 
high mass systems (M DD =2.4M ) with 0.7 < A/R Q < 4.2. In 
the next section the outcome of the close binary evolution is 
examined in connection to the possibility of producing final 
separations in the reference range 0.5 < A/R Q < 4.5, which 
corresponds to the relevant range of gravitational delays, for 
the full mass range of the DD systems. 

4.2. Shrinkage during the mass transfer phases 

The mass transfer phase in a close binary system may be dy- 
namically stable or unstable: in the first case the outcome is 
a wide system (occasionally wider than the primordial separa- 
tion), and the secondary may have accreted some of the donor's 
envelope mass. In the second case, a CE occurs and (generally) 
the system shrinks. The occurrence of one evolutionary chan- 
nel rather than another depends on the configuration of the ini- 
tial binary (e.g. the mass ratio at RLO, whether the envelope 
of the donor is radiative or convective, and so on). Anyway, 
the interesting systems are those which suffer a substantial de- 
gree of shrinkage during their evolution: the initial separations 
of the double CO WD progenitors range roughly from 100 to 
1000 R e , so that the first RLO takes place at all (upper limit), 
and that it does so after core Helium exhaustion (lower limit), 
when the star is on the Asymptotic Giant Branch. Thus, in or- 
der to merge within a Hubble time, the binary evolution should 



produce a total shrinkage on the order of ~ 10" 2 , 10~ 3 . This can 
be accomplished through one or more CE phases. 

The standard CE recipe (e.g. Webbink 1984) relates the ini- 
tial and final values of the binary parameters by requiring that 
the variation of the orbital energy is proportional to the binding 
energy of the envelope of the donor : 



( 24 ) m d (m d - m d ) 

1 v 1 f ' 



R 



m d m mf m 



2A f 



2Ai 



(26) 



where mf, mf are the mass of the donor respectively before and 
after the CE, m is the mass of the companion, A; and Af are 
the separations before and after the CE, R is the radius of the 
donor at contact (i.e. the Roche Lobe radius). a ce is parameter 
roughly describing the efficiency with which the orbital energy 
of the binary is used to expel the CE 5 : when smaller than unity 
the process is very inefficient, and the system emerges from the 
CE with a small separation; when a ce > 1 other energy sources, 
besides the orbital one, are used to expel the CE. Typical val- 
ues considered in the literature, and supported by hydrodynam- 
ical computations, range from Q- C e=0.5 to a ce =2 (Rasio & Livio 
[19961 . 

However, this formulation fails to explain the observed bi- 
nary parameters of 3 double He WDs, and Nelemans et al. 
(2000 1 propose an alternative equation for the first mass trans- 
fer, which consists in a parametrization of the system's angular 
momentum loss: 



Ji - J f AM 
— - — = YCE-rr- 

Ji M B 



(27) 



where 7i,7f are the initial and final angular momenta, AM is 
the mass lost from the system (roughly equal to the donor's 
envelope mass), Mb is the total mass of the binary before RLO, 
and ycE is a parameter, for which Nelemans et al. (|2000 1 find 
a value of ^ 1.5 by fitting the data. This formulation has been 
recently shown to describe the properties of a larger sample of 
binary WDs (Nelemans & Tout l2"00"5V 

Ea.J26> leads to a dramatic shrinkage, the typical ratio be- 
tween the separations after and before the CE phase (Af/A;) 
being on the order of a few 10~ 2 ; instead, for systems with a 
high mass ratio q (= m.2lm\), Eq. (I27> leads to modest shrink- 
age (Af/A; ~ 1). At low q, Eqs. i26\ and J27> yield similar 
values, and actually, Eq. (1271 is only applicable at relatively 
large q (see Appendix IA.2> . As for the second mass transfer 
episode, Eq. i26\ is generally adopted in the literature. I have 
then explored the results of two different evolutionary schemes: 
(i) the two successive RLO are regulated by Eq. i26\ : (ii) the 
first mass transfer occurs according to Eq. (I27> (when appli- 
cable, see Appendix IA.2l for a detailed description), while the 
second mass transfer is regulated by Eq. i26\ . The evolution- 
ary scheme (i) corresponds to the prescriptions used by e.g. 



Tutukov & Yungelson ( 1994 1, Ruiz-Lapuente & Canal ( 1995 



Han (1995), Han et al. (1998), Yungelson & Livio (2000i; the 
evolutionary scheme (ii) corresponds to the prescriptions used 
in Nelemans et al. (2001 1. 



5 Different formulations are found in the literature for Eq. 1261 . and 
correspondingly different meanings for a ce . I adopt here the formula- 



tion in Nelemans et al. < 



. with the geometrical parameter A=l. 
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Fig. 4. Exploration of the outcome of two mass transfer episodes in binaries with intermediate mass. The ratio between the 
separations of the final DD and the primordial systems is plotted as a function of the initial mass of the secondary for systems 
with mi > ni2- Eq. d!9l >. and the Eggleton ( 1998 1 formula for the Roche Lobe radius have been used. Blue and cyan points result 
from the application of Eq. J26t to both mass transfers; red and magenta points from applying either Eq. J27l > or Eq. J26t for the 
first mass transfer, depending on which of the two yields the wider separation At, At given m%, the larger m\ the smaller A/Aq 
is. Each panel is labeled with the adopted values for the parameters in Eqs. \26\ and J27i . The results in the top-left panel are 
reproduced in the other panels as small dots for an easy comparison. 



The total shrinkage, that is the ratio (A/Aq) between the 
separation of the newborn DD system and the original separa- 
tion of the binary, is shown in Fig.|4]as a function of the mass 
of the secondary component in the primordial system: blue and 
cyan dots refer to the (i), red and magenta dots to the (ii) evo- 
lutionary schemes. At fixed the points show the effect of 
varying m\\ the more massive m\ the smaller the A/Ao ratio. 
The four panels show the effect of varying the parameters a ce 



and ycE in Eqs. J26l > and ( 1271 . In order to gauge the effect of 
mass loss through stellar wind, which may occur prior to the 
CE phases, an additional parameter (/ en ) has been considered, 
such that the mass of the donor at each RLO is equal to a frac- 
tion / en of its initial mass. The sensitivity of the final shrinkage 
to this parameter proves significant. 

The two options lead to vastly different situations: for the 
scheme (i) the A/Ao ratio appears confined in a narrow range, 
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around a mean value which depends on a ce and / en . In addi- 
tion, there is a clear trend of A/Ao decreasing as mi increases. 
Indeed, more energy is required to expel the more massive CE 
in the more massive binaries. At fixed ar ce this implies that the 
systems shrink more and end up with a smaller A /Aq ratios. The 
(ii) prescriptions, instead, produce a wide range of A/Aq ratios 
at almost every m^. This means that systems with the same mi 
and Aq can end up very wide or very close depending on the 
mass of the companion. 

Turning now to consider the quantitative value of A/Aq, the 
three lines in Fig. |4]show the levels A/A = 5 • 1(T 4 , 5 • 1(T 3 , 5 • 
10~ 2 . Recalling that the initial separations of the double CO 
WD progenitors range roughly from 100 to 1000 R e , such lev- 
els correspond respectively to final separations of ^ 0.05,0.5 
and 5 R Q for the initially closest systems; to A 0.5,5 and 
50 R e for the initially widest ones. Inspection of Fig.|4]shows 
that the evolutionary scheme (ii) is capable of producing fi- 
nal separations in a very wide range, well including the range 
0.5R e $ A $ 5R e , leading to merging within a Hubble time. 
On the contrary, the scheme (i) appears to produce very small 
A/Ao ratios, so that only the initially widest binaries manage 
to merge on timescales on the order of some Gyr. In addition, 
the correlation between A/Ao and ma implies that lower mass 
systems have longer gravitational delays, as will be better illus- 
trated in the next section. 

It's important to notice the high sensitivity of the A/Ao ratio 
on the parameters (a ce ,f en ,ycE)', this, coupled with the high sen- 
sitivity of r gw on A suggests that the results of the binary evo- 
lution from the population synthesis codes are very dependent 
of the exact recipe used. At the same time, the correspondence 
between the {Aq, m\,mj) and (A, otir, W2r) is likely to be rather 
loose, and even more so if a distribution of (a ce ,/en,7cE) val- 
ues is realized in nature. In this respect, it is worth recalling that 
the computations here adopt a unique initial-final mass relation, 
while in reality, at a given initial mass, the remnant mass spans 
a (small) range, in relation to the precise point in the evolu- 
tion at which the mass transfer takes place. Therefore, for each 
pair (mi, mi) there will be a distribution of A/Ao, around the 
corresponding point in Fig. [4] 

This exploration of the results of the CE phases suggests 
to consider the two following extreme characterizations for the 
n(r gw ) function: 

- either, irrespectively of ;«2, the close binary evolution pro- 
duces a wide distribution of A and Mod, and these two vari- 
ables are virtually independent; 

- or the close binary evolution leads to a narrow distribution 
of the A /Aq ratio, so that the initially closest binaries merge 
in a short time, and the initially widest binaries tend to pop- 
ulate the long r gw tail of the distribution. In addition, the 
most massive binaries tend to end up with the smallest final 
separation, hence merge more quickly. 

The first characterization is more appropriate for the scheme 
(ii) of evolution, and will be referred to as WIDE DDs ; the sec- 
ond, suggested by the results of the evolutionary scheme (i), 
will be referred to as CLOSE DDs. Since both schemes pro- 
vide very small A/Aq ratios for some values of the parameters 
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Fig. 5. The g(r, r„) function for WIDE DDs, with f^ 2 =f^, for 
three values of the yS a parameter. The functions, plotted on arbi- 
trary units, are shown forr= 0.1, 0.3, 0.5, 0.7, 0.9, 1, 2, 5, 7 and 
10 Gyr. A very small minimum gravitational delay T gWi ; =0.001 
Gyr has been used. 

of the primordial binary system, the minimum gravitational de- 
lay r gw ; is likely to be very short. The dependence of the distri- 
bution function of the delay times on the parameter r gw i will be 
explored later, and only for the WIDE DD scheme. The maxi- 
mum gravitational delay r gw x is instead an important parame- 
ter for the characterization of the two schemes, as shown in the 
next sections. 



4.3. The distribution g(r, r n ) 

The fraction of systems which, having a nuclear timescale T n 
manage to merge within a timescale shorter than r gw = r - r n 
depends on the distribution of the gravitational delays. In this 
section the expressions for the g{r, r n ) function are derived for 
the WIDE DD and CLOSE DD schemes , based on the two 
different characterizations. 



4.3.1. WIDE DDs 

For the WIDE DD scheme the assumptions are: 

- Mud and A are independent variables; 

- the minimum gravitational delay r gw ; is independent of m%; 

- the maximum gravitational delay r gw x is larger than t - r n 
for all r n (at least for a total delay time up to the Hubble 
time); 

- the distribution function of the final separations follows a 
power law: n(A) oc A^\ 



Under these conditions, in Appendix IA.3l the following for- 
mula is derived: 



g(T, T n ) OC , w 





J\,2 



(T- 



11 ' gw,i 



for T< T\ 

for T\ <T< T2 



(28) 



where fi- d = 0.25(1 + y3 a ), T\ = r n + r gWji , r 2 = r n + r gWjX and 
is a function of r n , for which two cases are considered: 



'1,2 



/, W = 



VV0.75+0.75A, 
DD 



M, 



1.75+0.75/} a 
DD,x 



- M 



1.75+0.75j8„ 
DD,n 



(29) 
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with 

M D o,n = max(1.4,m 2R + 0.6); M DD;X = m 2R + 1.2; 
M DD = 1.4 + (m 2 -2)/6 

and W2R given by Eq. d!9i . The first case in Eq. &29\ ap- 
plies when there is a tight correlation between m.2 and its Mud 
progeny; the second applies when a wide distribution of Muu 
is obtained from systems with the same m 2 . 

Fig- |U illustrates the g(r, r n ) function (computed with the 
/* factor). In general, as r increases, the g function increases 
at every r n , since at longer delays more systems fulfill the con- 
dition r gw + r n < t at each t„. At the same time, for increas- 
ing t„ the number of systems with total delay smaller than a 
fixed t decreases, partly because of the decrease of the range 
of gravitational delays which fulfill the condition r gw < r — T n . 
A more thorough explanation of the trend of the g function can 
be found in Appendix IA.3I For j3 d = -0.9 most DDs are born 
with small separations, and the dependences on both r n and 
t are milder, since both the above mentioned effects are less 
relevant. Conversely, when /? a = +0.9 the distribution of the 
gravitational delays is skewed toward the large r gw values: the 
fraction of systems with total delay up to t greatly increases 
with t, as systems with longer r gw are included. 

It is worth pointing out that for f3 d < - 1 the g function is 
very sensitive to r gw ,j. This regime corresponds to distributions 
n(A) highly peaked at the minimum separation, a possibility 
that likely provides exceedingly small SNIa rates at late epochs 
in early type galaxies. 

4.3.2. CLOSE DDs 

The distinctive features of the end product of the CLOSE DD 
scheme of evolution appear to be that (i) the average A /Aq ra- 
tio is very small; (ii) it is correlated with the mass of the sec- 
ondary in the primordial system. As a result, the gravitational 
delays can be very small even for the systems with maximum 
primordial separation, especially for the most massive DDs. 

Fig. [6] shows the gravitational delay in this evolutionary 
scheme for systems born with Aq = 1 000 R & , and having adopted 
a ce =2 and / e n=0.9. Notice that these are about the widest inter- 
acting binaries: the components of systems born with a larger 
Ao evolve as single stars, and do not provide SNIa, at least 
in a Hubble time. In Fig. |6j for every m 2 (i.e. r n ), m t is de- 
creased from 8 to ra 2 6 , leading to longer gravitational delays 
(Mdd='«ir+'«2r decreases). The circled dots refer to systems 
with Mud ^ 1-4, and therefore suited to successful SNIa 
events. The upper envelope of the gravitational delays of the 
SNIa precursors is well represented by: 

log Tgw , x = min(- 16.66 + 3.17 log r n ; 6.02 + 0.52 logr n ) (30) 

with the delay times expressed in years. By varying the param- 
eters a ce and f en the linear regressions shift vertically, approx- 
imately maintaining their slopes. The vertical arrows drawn in 

6 Unlike for the WIDE DD scheme, it seems adequate to restrict 
here to systems with m lR > m 2 R. 




7.5 8 8.5 

lOg T n 



Fig. 6. Exploration of the gravitational delay of DDs from sys- 
tems with a primordial separation of Ao=1000 R Q , evolved 
through two CE phases regulated by Eq. (I26> . as a function 
of the MS lifetime of the secondary. Both ages are expressed 
in years. Circled dots mark systems with Mdd ^ 1.4M Q . The 
dashed line shows the 15 Gyr level; the arrows show the av- 
erage shift which applies when, leaving the other parameters 
unchanged, one adopts in turn different values for / en , a ce or 
Ao as specified. 

Fig-EDshow the amount of this shift for the indicated value of 
the parameters. 

The locus r gw = 15 Gyr is shown as a dashed line in Fig.|S] 
it appears that the possibility of realizing gravitational delays 
as long as the Hubble time is related to a fine tuning of the 
involved parameters, so that the shrinkage due to the two CE 
phases is not too severe. The systematic increase of the maxi- 
mum gravitational delay with increasing nuclear lifetime of the 
secondary reflects the smaller shrink of systems with smaller 
m 2 , which is related to the minor amount of energy required 
to expel a less massive CE. For the same reason, at given m 2 , 
a less massive primary implies a smaller amount of shrinkage 
at the first CE. Therefore, for the CLOSE DD scheme it seems 
appropriate to adopt a parametrization which emphasizes the 
systematic depletion of systems with long T ew as r n decreases, 
or the increase of the maximum gravitational delay in systems 
with smaller secondary mass. For simplicity, I consider directly 
r gw as the independent variable, and assume that for each m 2 
the differential distribution of the gravitational delays scales as 
n(r gw ) oc t™ between a minimum (r gw i) and a maximum value 
(i" gw ,x)- Admittedly, this choice is rather arbitrary; I just notice 
that, if the distribution of r gw is mainly controlled by the distri- 
bution of A, one can write: 

n(T gw )dr gw oc «(A)dA 
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Fig. 7. Illustration of how the t* point, defining the branching 
in the CLOSE DD distribution of gravitational delays, shifts 
with increasing total delay t. The straight lines show the loci 
y = t - r n for a few values of r; the curve shows the lo- 
cus T gWjX (r n ): for illustrative purposes, the relation log r gw x 
1 1 ii ti [ 16.66 + 3.17 log T n ; 5.4 + 0.52 log r n ] is used. Delay 
times are in Gyr. The asterisk marks the branching occurring 
at T n =r*. The shaded region shows the domain over which 
g(r,T n )=l. The dashed lines indicate the total range of nuclear 
delays, corresponding to the range of m.2 considered here, i.e. 
from 8 to 2 M@. 



which, in combination with the relations r gw oc A 4 and n(A) oc 
A A becomes n(T gw ) oc t^, with /J g = -0.75 + 0.25y6 a - 

With this assumption, the fraction of systems which man- 
age to merge within (r-T n ) is: 



ready exploded 1 . Thus, over the shaded portion of the plane, 
g(T,T n ) = 1, i.e. the fraction of systems which, having a nu- 
clear delay r n , merge within a total delay r is constant and 
equal to 1 . Clearly, at long total delays, the area in the parame- 
ter space over which systems have already merged widens, and 
only those with long t„, and therefore long r gw , still give a vari- 
able contribution. When r > T n x + T gw>x (T n , x ), g(r, T n ) = 1 for 
whatever r n , that is for all the SNIa precursors. Thus, in this 
scheme, the maximum delay time of the whole DD population 
from an instantaneous burst of SF is r x = T n , x +T gWiX (T n , x ). 

By varying the parameters of the close binary evolution 
(e.g. a ce , /en) the T gWiX (r n ) locus shifts. If evolution leads to 
systems which are all very close, the maximum gravitational 
delay will be small, and soon t x will be reached, at which point 
the distribution function of the delay times for the CLOSE DDs 
drops to 0. 

4.4. The Distribution Function of the Delay Times for 
DDs 

The distribution function of the delay times for the DDs is ob- 
tained by computing the derivative of Ea. J24> : 



r 1 



^ /-»min(T nx ,r) 



"(T n ) g(r, r n ) dr n . 



(32) 



where n(r„) is the distribution function of the nuclear delays of 
the SNIa progenitors, that is the function derived in Sect. 
3.2. 

Ea.P2> is solved by applying the Leibniz integral rule, as 
shown in Appendix I A. 41 The final function is: 



min(r llA ,r) 



«(t„) S w (t, t„) dr n for WIDE DDs (33) 



min(r n , x ,r) 



•JT n , 



complemented with 



«(Tn) S c (t, T n ) dr„ for CLOSE DDs (34) 



g(T, T n ) 



f o 



1 



' g».» ' gw.i 



for t - t„ < T gw>i 
for Tg W [ < t — T n ^ T gW)X 

for T - T n > T gw , x (T n ) 



(31) 



where T gw x is an increasing function of t„. It is further assumed 
that the minimum gravitational delay T gWj ; is independent of m% 
and small (compared to T n i), while the maximum gravitational 
delay is correlated to t„ via Eq. d30i . 

To better understand the behavior of the g function, Fig.0 
illustrates how the second branching in Eq. (13 1 1 varies as T n 
increases. The parallel straight lines show the y = t - t„ loci 
for selected values of t, ranging from 1 to ~ 18 Gyr; the curve 
is one example for the t„ w x locus (different from Eq. (13 0> for 
illustrative purposes). The intercept between each straight line 



and the T gw x curve (in T n =T*) defines the branching of Ea.l BTV 
at given t, systems with T n shorter than t* have all merged, 
while only a fraction of the systems with t„ > t* have al- 



/r(r)=o 

where 



5 W (T,T„) 



S c (t,t„) 



for t < t; and for t > t x 



7u(T-T„r a75+ °- 25A forr n <T-T w 

for T„ > T - T„ w ; 



.1+jSg i+ft 



for T n < T - T g 
for T n > T - T g 



Tn.inf — 



T n ,i for T < T n ,i + T gw , x (Tn,i) 



T* for T > TnJ + T gWjX (T n ,i) 

and with the following meaning of the symbols: 



(35) 



(36) 



(37) 



(38) 



7 Notice that for r < r n i + r gwx (T n i ) there's no intersection between 
the two loci: in this range, the total delay considered is so short that 
g(r, t„) < 1 even for the shortest r n . This regime is not visible in Fig.0 
because, for the chosen relation, T n>i + r gwx (T ni ) =0.07 Gyr. 
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- T"n,i, T n y. respectively the nuclear timescales of the most 
and least massive secondary in the SNIa progenitors' bi- 
nary systems; if ni2 ranges between 8 and 2 M Q , r n j 0.04 
Gyr and t„ x =* 1 Gyr; 

- r gWi i, T gw x : respectively the minimum and maximum gravi- 
tational delay of SNIa DD precursors, originated from sys- 
tems with given ;«2; r gw ,i is assumed independent of mi\ 
T gw ,x is assumed (i) larger than the Hubble time for all SNIa 
precursors in the WIDE DD scheme, (ii) correlated with r n 
as given by e.g. Ea. (l30> in the CLOSE DD scheme; 

- t* is the solution of the equation t - r n = T gw x (r n ), i.e. the 
minimum r n which contributes to the explosions at epoch 
r, which increases with r (see Fig.0; 

- Tj : the minimum total delay time (Tj=T n j+T gWj i), assumed 
to be independent of mi ; 

- t x : the maximum total delay time, larger than the Hubble 
time for the WIDE DD; equal to the maximum delay of 
the least massive SNIa progenitor for the CLOSE DD: 

Tx = Tn,x "t"Tgw,x (fn,x) > 

- : the term describing the dependence on the mass of the 
DD systems, given by Ea.( l29t ; 

- /3 a , j6 g : the exponents of the power law distributions 
adopted respectively for the final separations in the WIDE 
DD scheme, and for the gravitational delays in the CLOSE 
DD scheme; a flat distribution of A corresponds to /3 d = 0, 
and j6 g in the vicinity of -0.75 (if Mud varies in a small 
range). 

In addition, the (a, y) parameters need to be specified, in 
order to compute the distribution of the nuclear delays. In the 
models shown here, I have considered a = 2.35 and y = 1: the 
dependence of n(r n ) on these parameters is modest for nuclear 
delays up to ~ 1 Gyr (see Fig.[3j. 

In both Eqs @3 and at eacn tota l delay r, the f® D 
function results from the sum of the contributions from systems 
with a range of m%, and each contribution is proportional to a 
power of (t - r n ), the latter being just the gravitational delay 
of the progeny of systems born with secondary mass ma, which 
end up in a SNIa event at epoch r. 

The two equations correspond to different characterizations 
of SNIa precursors: Eq. 1331 separately accounts for the sen- 
sitivity of the gravitational delay from the total mass of the 
DD systems (through the factor ), and from the distribution 
function of the separations (through j8 a ). At the same time, it 
assumes that at any r up to the Hubble time, SNIa explosions 
come from all systems with r n < r n x . On the contrary, Eq. J34l > 
emphasizes the systematics of the range of the gravitational de- 
lays as 1112 (and hence r n ) varies: at each r, only systems with 
T n > r* contribute to the explosions, since those systems with 
shorter nuclear timescales have a too short maximum gravita- 
tional delay. 

Neither of the two relations will strictly apply in nature, but 
they can be used to investigate on the general shape of the dis- 
tribution function of the delay times in two extreme situations 
for what concerns the product of the Common Envelope evolu- 
tion. As an illustration, Fig.[8]shows the distribution function of 
the delay times for the DD model for one specific choice of the 
parameters, as labeled. For comparison, the distribution func- 
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Fig. 8. Illustration of the distribution function of the delay times 
in the DD model for the labeled values of the parameters. The 
red lines show the result for the CLOSE DD scheme; the blue 
and green for the WIDE DD scheme, respectively when Mod 
is strongly and mildly correlated with mo. Ea. P0> has been 
used for r gw x (T n ). For comparison, the magenta curves show 
the distribution function of the delay times in the SD model, 
for Chandrasekhar explosions. In all models the distribution of 
the nuclear delays is derived with a = 2.35 and y = 1. 

tion of the SD systems is also plotted. I recall here that /ia(i") is 
proportional to the SNIa rate following an instantaneous burst 
of SF. 

Similar to the case of the SD model, the distribution func- 
tion /j° d (t) appears characterized by three regimes: first a rapid 
increase, followed by a slow decrease or a wide maximum, 
and finally a late epoch, pronounced decline. This shape can 
be viewed as a modification of the function: for the DD 
model, early explosions are given by systems with short T n 
AND short r gw ; the flat portion corresponds to those epochs 
at which the SNIa events come from many combinations of 
r n and T gw ; at late epochs we are left with systems with long 
r gw . Notice that when r is large compared to T n>x Eq.(l33> can 
be approximated as /j° d (t) °c t - u - 75+0 - 25 A : compared to the 
SD model, this decline rate is considerably mild, and basically 
controlled by the dependence of the gravitational delay on the 
final separation A. A late epoch increase of f® D can be realized 
only if f3 d is large and positive, corresponding to n(A) distribu- 
tions highly skewed toward large values of the final separations, 
which is very unlikely. 

With respect to the WIDE DD, the CLOSE DD scheme 
yields a distribution of the delay times which is steeper both at 
the intermediate and at the late epochs. This reflects the relative 
paucity of systems with long r gw in this scheme of evolution. 
However, in spite of the very different assumptions, the overall 
behavior of the /j a functions in Fig. |S]look similar. In Sect. 5, 
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the difference between the models will be better quantified, by 
considering the distribution function of the delay times suitably 
normalized. 

For the DD model, at t= 1 Gyr there is a cusp: mathemati- 
cally this is due to the discontinuity of the S w and S c functions, 
coupled with the upper limit of integration for /j° D . In practice, 
the cusp occurs at the epoch at which the systems with smallest 
ni2 start contributing to the SNIa rate, that is T=T njX +T gWj i. After 
this epoch, increasing t corresponds to include systems with 
longer r gw but NOT longer r n . The prominence of the cusp is 
related to the f3 a ,(J3 g ) exponents which control the n(r ew ) distri- 
bution toward the minimum T gw i. 

Finally, the different scheme used to compute the distribu- 
tion function of m.2, i.e. whether using Ea. (ll6> or Ea. J20t has a 
very modest impact on the final ,/2? d (t) function. As for the SD 
model, the use of the Greggio & Renzini ( 1983 ) scheme yields 
a relatively larger rate at late epochs. 

4.5. Dependence on key parameters. 

I turn now to consider the dependence of the function f® D on 
the several parameters that need to be specified. 

Both the upper and the lower limits to the mass of the sec- 
ondary component in the SNIa progenitor systems are subject 
to some uncertainty, which reflects on the parameters T nj ; and 
r n x . The most massive CO DD systems come from progenitors 
in which both components are ~ 8 M Q stars; if these systems 
manage to produce a SNIa event, r n ,i is about 0.04 Gyr. As 
mentioned in the introduction, though, the ultimate fate of a 
double CO WD might be an accretion induced collapse, rather 
than a central Carbon deflagration, depending on the modali- 
ties of accretion (e.g. the accretion rate and the angular mo- 
mentum deposition on the WD). This question is highly de- 
bated in the current literature (e.g. Piersanti et al. 2003 and Saio 
and Nomoto 2005 1, to the aim of establishing the likelihood 
of the DD channel as SNIa precursors in general. However, 
notice that if the occurrence of the accretion induced collapse 
depends on the mass of the DD components, the DD channel 
remains a valid SNIa progenitor, but t^j changes to become 
the nuclear lifetime of the most massive secondary in a sys- 
tem which avoids the accretion induced collapse. Since in the 
current literature there's no claim of this effect, this possibil- 
ity is neglected here, and all the models adopt T n ;= 0.04 Gyr. 
As already mentioned, most double CO WDs come from sys- 
tems with 1112 greater than 2 M e , so that a reasonable value for 
r n x is 1 Gyr. However, the likelihood of a successful explo- 
sions may well be decreasing as mi approaches this limit, due 
to the requirement that Mod exceeds the Chandrasekhar mass. 
Therefore, r n x is treated as a parameter, and I show here the 
results obtained with r n>x =0.4,0.6 and 1 Gyr, corresponding to 
a lower limit to mi in SNIa progenitors of ~ 3, 2.5 and 2 M , 
respectively. 

Given the high degree of shrinkage which is obtained when 
applying the standard CE recipe, the lower limit to the gravita- 
tional delay (T gWj i) is likely to be very small. For this reason I 
adopt here a nominal value of r gWj j= 0.001 Gyr in most com- 
putations. However, the minimum gravitational delay could be 
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Fig. 9. Sensitivity of the distribution function of the delay times 
for DD progenitors on various parameters. The top panels refer 
to the CLOSE DD scheme, the bottom panels to the WIDE DD 
scheme of evolution. The left panels show the dependence on 
the T n , x parameter, with the color encoding as labeled. In addi- 
tion, the magenta line in the top panel has been obtained with 
T gw x (T n ) as given by Eq.j30> decreased by 0.6 Dex, mimicking 
the effect of lower a ce . The magenta line in the bottom panel 
shows instead the effect of adopting T gw ,i=0. 1 Gyr in the WIDE 
DD scheme. In all other models T gw ;=l Myr is used. The right 
panels show the effect of varying yS a and /3 g . The dotted lines in 
the bottom panels show the results of using =/*in the com- 
putation of /[° DW . All curves have been computed with n(r n ) 
derived from Eq.(I16> with (a, y) = (2.35, 1). 



larger, especially for systems with high initial mass ratio, if e.g. 
Nelemans et al. (2001 1 scheme of evolution applies. Thus, I ex- 
plore the sensitivity of the results on r gw ; only for the WIDE 
DD scheme, adopting a very large value of r gw i= 0.1 Gyr. For 
the CLOSE DD scheme the critical gravitational timescale is 
instead r gw , x : e.g. low values of a ce produce a maximum gravi- 
tational delay shorter than the Hubble time even for the lowest 
mass systems (see Fig. [6). Two options for the relation T gWiX (r n ) 
are considered. 

Without specific population synthesis computations, little 
can be said about the j3 d and the f3 g parameters. Since the dis- 
tribution of the separations in primordial binaries is typically 
taken to scale as n(A Q ) oc A Q l (e.g. Iben and Tutukov 1984, 
Han 1998, Nelemans et al. 2001), I consider the three values 
yS a = -0.9, 0, +0.9 which correspond to assuming that, as a re- 
sult of evolution, the distribution of the separations (i) remains 
basically unchanged; (ii) flattens off, so that any value of A is 
equally probable; (iii) changes slope, so that more DD systems 
are found with large A. The values of the j3 g parameter explored 
here are related to the three /? a values via fi g = -0.75 + 0.25y8 a 
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(see Sect. 4.3.2), and are f3 g = -0.975,0.75,-0.525. Notice 
that, since most combinations of (Mod, A) lead to short r gw , 
positive values for f3 g or even a flat distribution of the gravita- 
tional delays are extremely unlikely. 

Fig- illustrates the dependence of /^ d (t) on the param- 
eters related to the timescales (left panels), and on those re- 
lated to the distribution of the separations (right panels). The 
left panels show how the late epoch decline starts at r=T n x + 
T gw> i, as argued in the previous section. The more massive is the 
lower limit to 1112 of the SNIa progenitors' , the larger is the frac- 
tion of early explosions, i.e. the shorter will be the timescale for 
the release of the bulk of the nucleosynthetic products to the in- 
terstellar medium. If SNIas come mostly from DDs which are 
born with relatively wide separations, such that T gWj ; is large, 
the distribution function of the delay times behaves like the ma- 
genta curve in the lower left panel: the first explosion occurs at 
t = T"n,i + r gw after which the rate increases rapidly up to the 
start of the late epoch decline. The magenta curve in the top left 
panel shows instead what happens if the close binary evolution 
produces a large degree of shrinkage: the f® D function is more 
skewed toward short delays, and it does not provide systems 
with delays longer than r nx + T gWiX (T nx ). This illustrates the 
potential difficulty of accounting for SNIa in Elliptical galax- 
ies, if the close binary evolution produces too close DDs. 

The right panels in Fig. |9] show instead how the distribu- 
tion function of the delay times depends on yS a and on f3 g in the 
special case of T n , x =0.4 Gyr and r gw> i=l Myr. In general, f® D 
appears to be fairly sensitive to these exponents, with more and 
more early explosions as the distribution function of the sepa- 
rations of the DDs is more populated at the low A end. 

5. Discussion 

Very schematically, the distribution function of the delay times 
of SNIa progenitors derived in the previous sections for both 
Single and Double Degenerate models is characterized by: 

- an early steep rise; 

- an intermediate phase, hereafter referred to as the wide 
maximum; 

- a decline phase. 

The minimum and maximum delay times, the duration of the 
intermediate phase, and the slopes of the intermediate and late 
phases are different for the Single and the Double Degenerate 
progenitors, and are controlled by a few key parameters. 

For the SD progenitors, the minimum delay time is equal to 
the MS lifetime of the most massive secondary in the primor- 
dial binary producing a SNIa event; the wide maximum phase 
lasts until a delay time equal to the MS lifetime of the least mas- 
sive primary evolving into a CO WD suitable as SNIa precur- 
sor; the decline phase becomes very steep at late epochs, if the 
requirement of building up to the Chandrasekhar mass limits 
the SNIa progenitors to systems with more massive primaries 
in combination with less massive secondaries. In the decline 
phase, the slope of the fi a function depends on the IMF, and on 
the distribution of the mass ratios. In particular, the flatter the 
distribution of the mass ratios, the larger the fraction of systems 



with long delays, at fixed IMF slope. For the Chandrasekhar ex- 
ploders, the maximum delay is equal to the MS lifetime of the 
secondary whose envelope is massive enough to ensure that the 
most massive CO WD reaches the Chandrasekhar limit upon 
accretion. For the Sub-Chandra exploders, the maximum de- 
lay is equal to the MS lifetime of the secondary with a massive 
enough envelope to provide the minimum layer for Helium ig- 
nition on top of the companion. Both these constraints depend 
upon the efficiency of the accretion process. 

For the DD progenitors, the minimum delay time is equal to 
the MS lifetime of the most massive secondary in a SNIa pro- 
genitor system, plus the minimum gravitational delay; the wide 
maximum phase lasts up to a delay equal to the MS lifetime of 
the least massive secondary in a SNIa progenitor system, again 
plus the minimum gravitational delay. The slope of the decline 
phase is sensitive to the distribution function of the separation 
of the DD systems at birth. In addition, the overall distribution 
function of the delay times is steeper if a correlation exists such 
that the more massive binaries merge on a shorter timescale 
than the less massive ones, due to a more pronounced shrink- 
ing of the system. This happens in the standard treatment of 
the CE evolution (CLOSE DD scheme). Finally, the maximum 
delay time for the DDs is basically equal to the gravitational 
delay of the least massive and widest DD progenitor: if the CE 
stages were to induce a high degree of shrinking, the maximum 
delay could well be shorter than the Hubble time. 

5.1. Comparison with the results of populations 
synthesis codes 

These fi a functions have been derived with the aim of provid- 
ing a general characterization of the distribution function of 
the delay times for the various potential SNIa progenitors, and 
a number of convenient, though astrophysically motivated, ap- 
proximations have been introduced. It is thus very important 
to compare the general shape of these functions to the results 
of the population synthesis codes, which follow the individual 
evolution of close binaries in detail. Unfortunately, the distribu- 
tion of the delay times of the SNIa events, or equivalently the 
SNIa rate following an instantaneous burst of Star Formation 
(see Eq. 0), is not commonly found in the literature. Most au- 
thors rather quote the current SNIa rate in the Galaxy, which, 
following Eq. (|8j, gives information on the total realization 
probability of the la channel (i.e. Ai a ), but not on the shape of 
the /i a function. The most suitable paper to perform a detailed 
comparison between the analytic functions presented here and 
the results of a population synthesis code is the one by YL. 

Fig. 2 in YL shows the SNIa rate following an instanta- 
neous burst of Star Formation. Four types of precursors ap- 
pear in this figure: the DD-Ch, i.e. Chandrasekhar Double 
Degenerate exploders; the SG-Ch, produced by the evolution- 
ary path of the SD-Chandra considered here; the He-ELD and 
SG-ELD, which are two flavors of the Sub-Chandra Single 
Degenerate channel, the difference being that the former come 
from systems with 1112 greater than 2.5 M Q , which donate 
Helium to the degenerate companion, while the latter are sys- 
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terns in which mi is smaller than 2.5 M Q , which donate H, 
converted to Helium of top of the CO WD. The He-ELD and 
SG-ELD can be viewed as two complementary paths building 
up into the broad Sub-Chandra category considered here. 

The numerical simulations follow the evolution of one sin- 
gle population of binaries, which evolve through mutually ex- 
clusive channels; the fi a functions presented here, instead, are 
thought of as alternative to each other, for the total stellar popu- 
lation. One could consider a scenario in which SNIa come from 
different channels, and construct a composite analytic distribu- 
tion function of the delay times by properly assigning the var- 
ious key parameters, and the realization probabilities of each 
channel. I prefer to avoid this approach and perform the com- 
parison between the results presented here and YL's by taking 
into account the different mass ranges which evolve into the 
different channels. 

Reading off Fig. 2 in YL, the rate for the DD-Ch exploders 
starts at a delay time of Log t ^ 7.4, reaches a maximum shortly 
before 0. 1 Gyr, and then drops; for delays in excess of about 0.3 
Gyr the trend is close to a power law with a slope of ~ -1.2. 
At 10 Gyr, the rate has dropped of 2.2 Dex with respect to its 
value at maximum. YL state that in their simulations, the SNIas 
typically come from binaries with primary components in the 
range between 4 and 10 M , and that the DD-Ch channel ap- 
plies to the systems with secondaries more massive than 4 M . 
The upper limit on m\ is larger than the 8 M adopted here 
because the evolution in a close binary can prevent C ignition 
before the loss of the envelope in stars less massive than ~ 10 
M Q (Iben & Tutukov 1984). The evolutionary lifetime of a 10 
M star (with solar metallicity) is about 25 Myr, which is also 
the delay time at which the first DD-Ch events appear to oc- 
cur; the lifetime of a 4 M star is about 0.18 Gyr, close to the 
duration of the peak in the YL DD-Ch curve. Since YL adopt 
a description of the evolution during both CE phases similar to 
Eq. \26\ . the analogue of their DD-Ch case would be a CLOSE 
DD model with Tj= 0.025 Gyr and r n , x = 0.18 Gyr, while noth- 
ing can be said about the adequate value of the /3 g parameter. 
Using Eq. J30i with such short r n x , the maximum gravitational 
delay is much shorter than the Hubble time; on the other hand, 
as illustrated in Fig. |6j r gw x is very sensitive to the various 
parameters used to describe the CE evolution. In order to com- 
pare the analytic function for the DD model to YL results I 
consider a relation for r gw x obtained from Ea.J30> plus a zero 
point shift of -0.6, so as to recover maximum delays exceed- 
ing the Hubble time for the less massive SNIa progenitors. This 
case is shown in Fig.^|for three values of the j3„ parameter. It 
can be seen that the analytic fi a functions are very similar to the 
DD-Ch curves in YL; in particular, the case with f3 g =-0.75 is 
very well approximated by a power law with a slope of -1.27 
for t > 0.5 Gyr, and at 10 Gyr its /i a is 2.7 Dex lower than 
its maximum value. This is a remarkable similarity, given the 
completely different ways in which the two functions have been 
obtained. Notice that in the range 8.7 < log t < 9.5 the match 
is better than this, since the analytic curves steepen when ap- 
proaching the maximum delay time. 



To be precise, YL quote a secondary mass below (2-3 M Q ) 
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Fig. 10. Distribution function of the delay times for DD pro- 
genitors for a choice of parameters which represent the DD-Ch 
evolutionary channel in YL. See text for more details. 

It is worth to point out that the f\. d functions for the CLOSE 
DD scenario are also in broad agreement with the SNIa rate fol- 
lowing an instantaneous burst of star formation for other ren- 
ditions of the DD scenario (Tutukov & Yungelson 1994 Ruiz- 
Lapuente & Canal |T9*98 1. 

The similarity between the analytic and numerical func- 
tions for the SD, Sub-Chandra cases is much less apparent; still 
the differences can be understood in terms of different mass 
ranges populating this channel. After a steep rise, starting at 
log t = 7.5 (close to the MS lifetime of an 8 M star), the rate 
from the He-ELD channel has a wide maximum, followed by 
a rather abrupt drop, which sets in at about 0.5 Gyr. The lat- 
ter is the evolutionary lifetime of a ~ 2.7 M star, not far from 
the least massive secondary populating this channel, i.e. 2.5 
M . Systems with lower secondary mass do evolve on a longer 
timescale, but do not go through the He-ELD channel. On the 
other hand, the rate from the SG-ELD exploders starts at a de- 
lay time of 0.6 Gyr, in correspondence to the lifetime of the 
most massive secondary which evolves through this channel. 
The rate rapidly reaches a maximum and then starts declining 
with a trend close to a power law with a slope of -1.6. At a de- 
lay time of ~ 5.6 Gyr the rate drops dramatically. For compari- 
son, the slope after the wide maximum of the blue, dot-dashed 
y^SD curve m Figgis =i -1.7. It is tantalizing to conclude that 
the analytic formulation presented here for the Sub-Chandra 
exploders reproduces the general trend of the He-ELD + SG- 
ELD channels in YL, with the late dramatic drop possibly re- 
lated to inefficient accretion from low mass secondaries. 

Finally, the rate for the SG-Ch channel is indeed very 
different from the analogue analytic /? D functions in Fig. 
However, as stated by YL, only systems with mi lower than 
2.5 M evolve along this path, the more massive secondaries 
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becoming either a DD-Ch, or a He-ELD. Therefore, the late 
start of the rate for this channel is understood in terms of a late 
T n ;; the rapid drop at delays in excess of 1 Gyr could instead 
reflect a low accretion efficiency. Actually, the comparison of 
the predictions for the Single Degenerate Chandrasekhar chan- 
nel is particularly difficult because of the different approaches, 
as anticipated above: the systems evolving through the SG-Ch 
channel are clearly a minority in the population synthesis code. 

The general conclusion is that the analytic distribution 
function of the delay times presented here can provide an ex- 
cellent match to the results of the population synthesis codes, 
once the appropriate mass ranges and evolutionary timescales 
are assumed for the individual evolutionary paths. Obviously, 
the numerical simulations also estimate the realization proba- 
bility of the various channels, and the total realization probabil- 
ity of the la event, that is the Ai a factor. There's no attempt here 
at evaluating this factor, which can either be taken from obser- 
vational estimates, following Eq. (|8), or from the population 
synthesis results. 

The analytic approach offers several advantages, most no- 
tably: 

(i) an easy way to explore the consequences on the evolution of 
stellar systems of the various candidates, viewed as alternative 
SNIa progenitors; 

(ii) it has a built in parametrization of the key properties of the 
alternative candidates, i.e. mass ranges, IMF and distribution 
of the mass ratios, distribution of the separations of the DD 
systems; 

(iii) a flexible tool to build up an overall distribution of the 
delay times by mixing the different analytic /i a functions, 
weighted by the relative contribution of the individual channels 
to the total realization probability, as astrophysical considera- 
tions may suggest. 

5.2. Comparison between different SNIa candidates 

For a meaningful comparison of the distribution functions of 
the delay times from the various potential SNIa progenitors it 
is necessary to normalize the /j a functions. Among the var- 
ious possibilities, a convenient normalization is to consider 
X / /ia( T ) • dr = 1 (see Sect. 2): in this way, the specific SNIa 
rate at the current epoch in a system which formed its stars in 
an initial star formation episode of duration Af is given by (see 
Eq. OS): 



= K ■ Ala ■ </la>f-£ 



(39) 



M 

where (fia)t-&t,t is the average of the normalized distribution 
function of the delay times over the age range of the stars in the 
system, and Ai a is the realization probability of the SNIa sce- 
nario. I recall here that k a depends on the IMF, and is equal to 
1.55, 2.83 respectively for Rroupa and Salpeter IMFs. Fig. II II 
shows some selected cases of the analytic distribution functions 
normalized to 1 in the range between and 12 Gyr. Although 
the maximum delay for the different candidates plotted here is 
larger than 12 Gyr, most distributions steepen enough at late 
times that the fraction of SNIa explosions at delay times ex- 
ceeding 12 Gyr is negligible. In the attempt of showing a wide 




Fig. 11. Distribution function of the delay times for SD and 
for DD progenitors as labeled. The plotted functions have been 
normalized to 1 in the range r < 12 Gyr, and the normalization 
constant is in units of Gyr 1 . For the all models the distribution 
function of the secondaries has been derived using Eq.jl6> with 
a = 2.35. See text for more details. 



range of possible solutions, the plotted cases adopt the follow- 
ing parameters: 

(i) for the SD models, the solid lines have been computed with 
y=l, the dashed lines with y=0.005; the latter choice enhances 
the relative number of low mass secondaries (see Fig. [5}; 

(ii) for the DD models, the solid lines have been computed with 
#,=0, j6 g =-0.75; the dashed lines with jS a =+0.9,y6g=-0.975 for 
the WIDE DDs and the CLOSE DDs respectively. The latter 
cases are meant to illustrate the wide range of decline rates 
allowed by these options. 

For a given realization probability of the SNIa event (Ai a ) 
all models plotted in Fig. [^provide the same total number of 
SNIa within 12 Gyr from a burst of SF of given total mass, but 
the events are differently distributed in time. Although for all 
of the models most of the events occur within the first 1 Gyr, in 
the WIDE DD scheme the explosions are more evenly spread 
over the whole 12 Gyr range. Both the CLOSE DDs and the 
SD Sub-Chandra model yield an age distribution of the events 
very skewed at the early epochs; the SD-Chandra model ex- 
hibits a dramatic drop at late epochs, related to the requirement 
of building up to the Chandrasekhar mass by accreting the en- 
velope of low mass secondaries. 

The models in Fig. [^ can be compared to the observed 
SNIa rate per unit mass in elliptical galaxies measured by 
Mannucci et al. ( 120051 hereafter M2005), which is quoted of 
0.044 (+0.016)(-0.014) SNuM, or 0.044(+0.016)(-0.014) x 
1CT 3 events per M per Gyr. If 10% of the stars in the mass 
range from 3 to 8 M Q end up as SNIa, the factor {k a Ai a ) is 
about (2, 3) X 10~ 3 respectively for Salpeter and Kroupa IMF. 
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Adopting k a A\ a = 2.5 x 10 3 , and inserting (nf a (t)/M) = 
0.044 10~ 3 , Eq. (39} yields (fu) t - h „= 0.0176 Gyr" 1 , which is 
the level indicated by the black line in Fig. ^2 The dashed re- 
gion shows the upper and lower limit relative to the range of the 
SNIa rate quoted by M2005 for E/S0 galaxies. For a larger re- 
alization probability of the SNIa event, following either from a 
wider range of progenitor masses, or from a larger probability 
of the SNIa channel within a given mass range, the observa- 
tional constraint shifts downward. 

The intercept between the theoretical fi d functions and the 
observational constraint yields the average age which the stars 
in early type galaxies should have in order to reproduce the 
data. So, e.g. the SD-Chandra model is able to fit the data if 
either the stellar population in Es is young, or if the realization 
probability of the SNIa scenario is larger than what adopted. 
Both the age of the stars in Es and the k a Ai d factors are uncer- 
tain; therefore Fig.^Jdoes not allow us to draw stringent con- 
clusions about the best model for SNIa precursors. However, 
the figure shows the interplay between the various quantities. 

At 12 Gyr the SD models with y- \ fall short by about one 
order of magnitude with respect to the level indicated by the ob- 
servations. The mismatch is more severe for the Chandra case. 
Such a big discrepancy is difficult to recover either by increas- 
ing the realization probability of the SNIa scenario, or by in- 
voking a younger age for Ellipticals: on the one hand the black 
level in Fig. ^2 already assumes that an important fraction of 
stars (i.e. ^ 10 %) in the suitable mass range end up as SNIa. 
On the other hand, M2005 data refers to a sample of more than 
2000 early type objects, and the spectrophotometric properties 
of this class of galaxies strongly suggest that they are old (see 
e.g. Renzini 1999; Peebles 2002 1. The only possibility to rec- 
oncile the SD model with the SNIa rate in ellipticals seems to 
be that of assuming a very low y (thin curves), so that the dis- 
tribution of the secondaries is maximally populated at the low 
mass end. Even so, an accretion efficiency close to 100 % is 
required for the SD-Chandra models to meet the observations, 
which seems unlikely, as argued in Sect. 3.2. 

The DD models more easily account for the observations, 
provided that gravitational delays as long as the Hubble time 
are realized, i.e. the common envelope phases do not lead to a 
too severe shrinking of the DD systems. The steepest f\. d func- 
tion {CLOSE DD with r n , x =0.4 Gyr and f3 g =-0.975) fall short 
by a factor of ~ 5 with respect to the observational limit, and 
are thus unfavored. Notice that the mild slope of the fi a func- 
tion from intermediate ages onward, implies that assuming a 
younger age for Ellipticals does not efficiently improve the 
fit for this kind of models. The illustration clearly shows that 
lower Ai a and/or older ages for ellipticals are accommodated 
with WIDE DD models. 

The normalization chosen for Fig. ^2 corresponds to as- 
suming that all models, with the same realization probability, 
yield the same total number of SNIa out of a stellar generation 
of unit mass (and therefore the same chemical enrichment). In 
a different approach, Fig. ^] shows the models normalized to 
their value at 10 Gyr: this corresponds to forcing all of them 
to fit the current SNIa rate in Es with an average stellar age 
of 10 Gyr, albeit with different values for the factor £„Ai a . The 
realization probabilities required by this normalization are of 
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Fig. 12. Distribution functions of the delay times for the same 
models shown in Fig. ^2 with the same color and line coding. 
The functions are here normalized to give the same value at 10 
Gyr in order to illustrate the different early epoch behavior of 
models which reproduce the current SNIa rate in Ellipticals. 

the order of 10" 2 for the SD models, of 10" 3 for the DD mod- 
els, but there is a noticeable dependence of the A Ia factors on 
the various parameters defining the models, including the for- 
malism to describe the binary population (i.e. whether Eq. d!6i 
or Eq. J20i are used). Only the SD-Chandra models with low 
e j do require a totally unrealistic realization probability, cor- 
responding to ^ 100% of the stars with mass between 2 and 
8 M Q , and will not be considered further. Once normalized in 
this way, the various models correspond to dramatically differ- 
ent evolution over cosmic time of the SNIa rate from a burst 
of SF (notice that the rate is plotted on a logarithmic scale). 
This property offers an important tool to discriminate among 
the models by looking at the impact on the large scales, like the 
Iron Mass-to-Light ratio in Clusters of galaxies, the evolution 
with redshift of the SNIa rate in Ellipticals, and the systematic 
trend of the SNIa rate with galaxy type (Greggio 2005 1. The 
latter point is addressed in the following section. 

To conclude this section, Fig. ^] shows the cumulative 
number of SNIa explosions as a function of age for an instan- 
taneous burst of SF, normalized to the total number of events 
within 12 Gyr. This figure illustrates another distinctive charac- 
teristic of the various SNIa models: the time scale over which 
50 % of the total SNIa explosions from an instantaneous burst 
of SF have occurred. This timescale can be taken as indicative 
of the typical delay with which the Fe from SNIa is released 
to the Interstellar Medium, thereby decreasing its a/Fe ratio. 
Within the explored range of models such timescales varies be- 
tween 0.3 and 3 Gyr going from the steepest CLOSE DD case, 
to the flattest WIDE DD case. The SD-Chandra model with 
7=0.005 (green dashed line) has a typical timescale of ^ 1 Gyr, 
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Fig. 13. Cumulative number of SNIa explosions following an 
instantaneous burst of SF for the same models shown in Fig. II II 
with the same color and line coding. In particular, the thin blue 
and cyan lines are especially flat WIDE DD cases; the thin red 
and magenta lines are especially steep CLOSE DD models. 



which is often taken as a reference value to infer the forma- 
tion timescale of systems exhibiting an a/Fe overabundance. 
Fig. 1131 emphasizes that such timescale does depends on the 
SNIa model, namely it is longer the flatter the distribution of 
the delay times of the SNIa progenitors is. The estimated for- 
mation timescales, then, remain uncertain by a factor of a few, 
modulo the actual SNIa channel that dominates in nature. 

5.3. SNIa rate in different galaxy types 

As anticipated in the previous section, the large differences of 
the SNIa rate temporal behavior typical of the various mod- 
els translate into a different trend of the SNIa rate as a func- 
tion of the galaxy type, thereby offering a tool to discriminate 
among the potential progenitors. This has already been outlined 
in Sect. 2, where the Cappellaro et al. ( 1999 1 data have been 
shown to indicate that the ratio between the /j a value at late de- 
lay times and its average value over the whole range of delay 
times (up to the Hubble time) should be 0.15. For the ana- 
lytic functions, the quantity is equal to 0.02 (0.03), 0.05 
(0.08) respectively for the SD Chandra (Sub-Chandra) models 
with y = 1 and 0.005 ; equal to 0.07 and 0.15 for the CLOSE 
DDs with T n>x =0.4,l (and y8 g =- 0.75); while it is 0.23 for the 
WIDE DDs with a flat distribution of the separations (B d =0). 
Therefore, the ratio between the SNIa rate in Ellipticals and 
Spirals indicates that the Single Degenerate model underesti- 
mates the current rate in early type galaxies, with respect to the 
rate in late types, a result of its fast decline at late times. 



This constraint has been derived considering a schematic 
description of the star formation history in early and late type 
galaxies, and by using a theoretical value for the M/Lb ra- 
tio in the two galaxy types. A much better constraint on the 
SNIa model progenitors can be built upon the recent results 
by M2005, by considering the trend of the SNIa rate per unit 
galaxy mass with the parent galaxy type. In fact, Eq. (|4) can be 
written as: 

T j^-=k a -A Ia -(f Ia ) m (40) 

where (/i a )^(«) is the average of the distribution function of the 
delay times weighted by the star formation rate over the galaxy 
lifetime, and where the integral of the SFR has been approx- 
imated with the galaxy (stellar mass) M 9 . Eq. J40i clearly 
shows that the trend of the SNIa rate with galaxy type reflects 
the systematics of the their star formation histories. Given that 
/i a is a decreasing function of the delay time, the younger sys- 
tems will have a higher SNIa rate per unit mass, by an amount 
which depends on the shape of the /i a function. The correlation 
in M2005 can be interpreted as the result of Eq. J40i . where the 
different SF histories imply both a different SNIa rate per unit 
mass, and a different B - K color. 

Figures ITU and [TBI show the M2005 observed correlation, 
and the theoretical predictions for the various SNIa models, 
computed as follows. To describe the SF history in the various 
galaxy types I have considered for four families of models: 

(1) old burst models, characterized by a constant SFR within a 
starting epoch t — and ending epoch t-g, = 0.5, 1, 1.5 and 2 
Gyr; 

(2) exponentially decreasing models, all starting at epoch t = 
and with different e-folding times ?sf: <A(0 = e "sf with fsF = 
1,3,5,7,9, 11 Gyr; 

(3) exponentially increasing models: all starting at epoch t — 

12-f 

and with a different characteristic time r+: tfr(t) — e r * with 
t* = 8, 6, 4, 3, 2, 1 Gyr. These models correspond to a sequence 
of age distributions all peaked at 12 Gyr, and decreasing e- 
folding age, so that the lower t* the younger the average age 
of the stellar population; 

(4) young burst models, characterized by a constant SFR within 
a starting epoch = 9.5, 10, 10.5, 11 Gyr and all ending at 
t = 12 Gyr. 

In addition, a model with a constant SFR over the whole range 
of 12 Gyr has been computed. The considered range of SF 
histories appears to encompass the range of B - K colors of 
the galaxies in M2005 sample. The composite colors have 
been computed using Girardi et al. (2000) simple stellar pop- 
ulation models with solar metallicity. The bluest and reddest 
bin in M2005 are defined as lower and upper limits only (i.e. 
B - K < 2.6 and B - K > 4.1); here they have been specified 
so that they are close to the bluest and the reddest model stellar 
population. The point type in Figs 1141 and 1 151 encodes the SF 
history: burst models are plotted as circles; exponential models 
as crosses, each family (2 and 3 described above) connected by 

9 As pointed out in Sect. 2, the actual stellar mass in a galaxy can 
be « 30 (45) % lower than the integrated SFR for Salpeter (Kroupa) 
IMF, due to the mass return form dying stars. 
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Fig. 14. Comparison between Double Degenerate model predictions and observations of the SNIa rate per unit mass as a function 
of the color of the parent galaxy. The data from M2005 paper are plotted as black dots, with their quoted error bars, and connected 
with a solid line. The 9 panels show the trend with B - K as a tracer of SF history (see text) of CLOSE DD (blue) and WIDE 
DD (red) models: from left to right the minimum secondary mass in SNIa progenitor systems decreases, and its MS lifetime 
T n ,x increases; from top to bottom the slope in the decline phase varies as indicated in the central panels, reflecting different 
distribution functions of the DD separations. Also labeled are the realization probabilities required by the different models in 
order to reproduce the current rate in the reddest galaxies. All models adopt a = 2.35 and y = 1, and both options for the 
distribution function of the secondaries, i.e. Ea.(ll6> and Ea. (l20> are displayed. The minimum gravitational delay is assumed to 
be negligible; the maximum gravitational delay for the CLOSE DDs follows Eq. ( I30> . The two options for the WIDE DD case, 
i.e. =/j W and /*in Eq. d36l give the same results, and indeed are both plotted in of the various panels. 



a line. The solitary dot shows the model with a constant SF rate 
over 12 Gyr. 

The theoretical rates have been normalized so as to repro- 
duce the observed values in the reddest galaxies with a SFH 
given by the old burst model with a duration of 2 Gyr; con- 



sequently, the normalization implies different realization prob- 
abilities of the SNIa event for the different models. The val- 
ues for Ai a required by this normalization, having adopted a 
Kroupa IMF and Eq. d!6i for the distribution function of the 
secondaries, are labeled in the figure. Ai a is the number frac- 
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Fig. 15. Comparison between Single Degenerate model pre- 
dictions and M2005. The point type encodes the SFH, as in 
Fig- El All models adopt a = 2.35; other parameters are la- 
beled. The blue and red curves refer to models computed with 
n(m2) given by Ea.(ll6>: the green and magenta loci show the 
results obtained with Ea. (l20> and y = 1. 

tion of SNIa events out of a stellar generation, and should be 
compared to the number fraction of stars in a mass range which 
can lead to the event, e.g. 3 < m/M Q < 8. For Kroupa IMF, 20 
% of the stars fall in this mass range, and therefore the normal- 
ization requires that in the WIDE DD scenario, approximately 3 
% of the stars in the suitable mass range should end up as SNIa. 
For the CLOSE DD scenario the analogous fraction is around 
~ 10 % . The dependence of these figures on the parameters for 
the DD model can be appreciated from the values labeled in 
Fig. 1141 for Salpeter IMF, the required fraction of SNIa events 
from the same mass range is smaller by a factor of ~ 0.6. It 
can be noticed that choosing any other old burst model, or the 
oldest exponentially decreasing model, would not change the 
normalization. 

It appears that the both families of DD models do fit well 
the observational data when the choice of the parameters is 
such to provide an intermediate shape of the distribution func- 
tion of the delay times. The CLOSE DDs yield a too steep evo- 
lution of the SNIa rate per unit mass with galaxy color, for a 
low jfij, and a short r n x ; similarly, the WIDE DDs give a too flat 
relation if f3 a is large, as well as t„ iX . This comparison does not 
necessarily favor the WIDE or the CLOSE DDs; rather it points 
to a moderate solution: either relatively flat CLOSE DDs or 
relatively steep WIDE DDs. 

Fig. ll5l shows the analogous results for the SD models. The 
left panel clearly shows that, with the standard choice of the 
parameters (i.e. n{m.2) as in the population synthesis computa- 
tions, a = 2.35 and y — 1) both Chandra and Sub-Chandra 
exploders imply a too large increase of the SNIa rate per unit 
mass going from early to late type galaxies. The right panel 
shows the results for a choice of the parameters aimed at im- 
proving the match between the SD model and the observations. 
When using Eq. (I20> to describe n(ni2) the agreement is bet- 
ter (green and magenta points), but is seems that only with 
a very low value of y (i.e. if all the mass ratios q = m2/m\ 
are equally probable) can the SD model be reconciled with 
the observations. In this case, the normalization to the rate in 
the reddest galaxies implies that ~ 10% of all stars born with 



2 < m/M Q < 8 should end up as SNIa of the SD Chandra vari- 
ety; alternatively ~ 12% of all stars born with 3 < m/M Q < 8 
should end up as SD Sub-Chandra SNIa. 

The comparison between the M2005 data and the theoret- 
ical models does not definitely rule out any of the alternative 
progenitors, but puts constraints on the key parameters within 
the various families. The precise values of these parameters 
are subject to some uncertainty. M2005 convert their observed 
SNIa rate per unit luminosity into a rate per unit mass by using 
the results from the galaxy models by Bell & de Jong (2001 1, 
which are much more complicated than the SF histories con- 
sidered here. This introduces a (small) inconsistency between 
the models and the data. Another caveat concerns the approx- 
imation of the galaxy mass with the integrated SF rate, which 
affects the zero point (by not more than a factor of 2), and the 
slope of the theoretical trend. Both these approximations, how- 
ever, hardly affect the major conclusion that indicates that 

(i) the distribution function of the delay times should drop by a 
factor of « 100 from its early peak to its value at 10 Gyr; 

(ii) in principle all candidates can reproduce the observed trend 
of the SNIa rate per unit mass with galaxy color, with an ad- 
equate combination of the key parameters. However, the SD 
channel seem to require quite some fine tuning; in particular 
the Chandra model is viable only if the accretion efficiency is 
very large, which appears unlikely. 

6. Summary and concluding remarks 

In this paper I have presented a straightforward formalism to 
relate the rate of SNIa events in stellar systems to their star for- 
mation history through two fundamental characteristic of the 
SNIa progenitor model: the realization probability of the SNIa 
scenario from a single age stellar population Ai a , and the distri- 
bution function of the delay times /i a (r), which is proportional 
to the SNIa rate past an instantaneous burst of star formation. 

The various models for the SNIa progenitors correspond to 
different values for A\. d , and to different shapes of the fi a func- 
tion. No attempt is made here to give a theoretical value for 
Ai a , which can be derived either from the numerical realiza- 
tions in population synthesis models, or directly from the ob- 
servations. The latter suggest Ai a on the order of 10~ 3 (see Sect. 
2.1). A more precise estimate could be derived with a detailed 
modeling of the star formation in late type galaxies, and with 
a robust assessment of the IMF. For the /i a function instead, 
the paper presents analytical formulations to describe the most 
popular SNIa progenitor models. It is shown that some param- 
eters play a key role in shaping the distribution function of the 
delay times, most notably: the mass range of the secondaries in 
systems which provide SNIa events; the minimum mass of the 
primary which yields a massive enough CO WD to ensure the 
explosion, and the efficiency of accretion on top of it, for the 
SD model; the distribution of the separations of the DD sys- 
tems at their birth. In addition, it is found that the algorithm 
adopted to describe the distribution of the stellar masses in bi- 
naries has an important impact on the slope of the /i a function 
for the SD model. 

In all of the three scenarios (SD Chandra and Sub-Chandra, 
and DD Chandra) the distribution function of the delay times is 
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characterized by a steep, early rise, so that the maximum rate is 
reached soon after the first SNIa explodes. For the SD model, 
this reflects the behavior of the clock (m^) which describes how 
the rate of change of the evolutionary mass decreases as the 
delay time increases. The quantity \niT\ drops so fast to pre- 
vail over the increase of the number of SNIa precursors as the 
delay time increases. For the DD model, the clock has an addi- 
tional contribution from the gravitational delay; still the early 
steep rise is present because most DD systems have a short r gw , 
given that a wide range of masses and separations correspond 
to a short gravitational delay. Therefore, the property of the 
distribution function of the delay times to reach its maximum 
shortly after the first event occurs has a very robust justification. 
As repeatedly said, following an instantaneous burst of SF, the 
epoch of the first event is equal to the MS lifetime of the most 
massive secondary (r n i), for the SD model, plus the minimum 
gravitational delay (r gw ;), for the DD model. There's no appar- 
ent reason to consider an upper limit for ;«2 in SNIa progenitor 
systems smaller than the most massive primary which provides 
a CO WD (i.e. an ~ 8 M star). This implies T n ,; 0.04 Gyr. 
As shown in Fig.|4] even in the WIDE DD scheme the binary 
WDs can emerge from the CE with a very small separation, 
so that also r gWi ; is likely to be very short. Therefore, in the 
context of the binary evolution, it seems very hard to accom- 
modate a shape for the distribution of the delay times similar to 
the one that according to Strolger et al. (2004) best explain the 
redshift dependence of the SNIa rates measured in GOODs, i.e. 
a Gaussian distribution centered on t ~ 3 - 4 Gyr. Such a dis- 
tribution could be obtained if SNIa events were produced only 
by systems with mi in the vicinity of 1 .3 M Q , for the SD model, 
or by DD systems born with a separation of about 2.5 R e (but 
notice that a little spread of the binary mass implies anyway 
a large spread of delays). Both hypothesis are very contrived. 
Furthermore, a Gaussian distribution would not meet the con- 
straint on the /i a function derived from the trend of the SNIa 
rate with parent galaxy type, like in M2005. 

The analytical functions derived in this paper compare very 
well to the results of the population synthesis codes in the lit- 
erature, and especially to those in Yungelson and Livio (2000), 
when taking into account the appropriate mass ranges and kind 
of progenitor. The Monte Carlo simulations follow the evolu- 
tion of a population of binaries, which evolve through mutually 
exclusive channels, according to the response of the system to 
the mass exchange phases. The remarkable similarity between 
the analytic functions derived here and the numerical results 
suggests that the shape of fi a is mainly determined by the mass 
ranges and general characteristics of the clock of the explosion, 
while the details of the response of the individual systems to the 
RLO events are of lesser importance. 

Once put on the same formalism, the various models for 
the SNIa progenitors can be compared to the relevant observa- 
tions in order to judge which one best accounts for the data. 
The most direct observational counterpart of the functions 
is the redshift dependence of the SNIa rate in elliptical galaxies 
(per unit galaxy mass), due to the fact that the bulk of their stars 
formed at very early epoch, and the single burst is a reasonable 
approximation. Along this line, Maoz & Gal- Yam ( 2004 1 have 
compared the SNIa rate in galaxy clusters at redshift between 



and 1 to the family of fin functions proposed by Madau et al. 
( 1998), reaching the conclusion that the average delay time of 
SNIa precursors ought to be < 2 Gyr. The notion of a typical 
delay time for the SNIa precursors has very little justification in 
the context of stellar evolution in binaries, since a wide range of 
delay times is produced by any kind of progenitor. Rather, the 
slow increase of the SNIa rate with redshift reported in Maoz 
& Gal- Yam ( 2004) would favor the DD model, perhaps better 
of the WIDE DD variety. With the progress of the Cluster SN 
surveys (see e.g. Maoz ( 2005 )) we will be able to further inves- 
tigate the redshift dependence of the rate of SNIa in Elliptical 
galaxies, and come to stronger conclusions. A less direct, but 
definitely complementary and effective approach to constrain 
to SNIa progenitor model is attempted here, by considering the 
predicted trend of the SNIa rate (per unit mass) with the parent 
galaxy type (see also Delia Valle & Livio 1994, Ruiz-Lapuente, 
Burkert & Canal [T995V The M2005 data suggest that the DD 
channel is favored with respect to the SD channel, and that 
the distribution of the separations of the DD systems should 
be such to produce a moderate decline of the fi a function at 
late times. The SD model is not completely ruled out by this 
comparison, but it requires (i) a flat distribution of the mass 
ratios and (ii) accretion efficiencies close to 100 %. This can 
be accomplished only if the matter is accreted and burned on 
top of the WD at the same pace, so as to avoid either expan- 
sion beyond the Roche Lobe (and the formation of a Common 
Envelope), or the accumulation of a Hydrogen layer which is 
eventually ignited under degenerate conditions (so that a Nova 
explosion occurs). Even in the Hachisu et al. ( 1996) scenario 
part of the accreted matter is lost by the system in a stellar 
wind. 

The various models for the SNIa progenitors have different 
impact on the large scales; some preliminary considerations are 
in Greggio ( 2005 ), while more detailed investigations will be 
presented elsewhere. Here I just remark a few points. 

Once normalized to reproduce the current SNIa rate in 
Ellipticals, the SD model corresponds to a large Fe mass to 
light ratio in Cluster of Galaxies (see Fig. 13, and Greggio 
2005). A detailed study of the expected Fe mass to light ra- 
tio in galaxy clusters as a function of the various possible pro- 
genitors, and its evolution with redshift will allow us to better 
constrain the SNIa model and the contribution of SNIa to the 
Fe in the intracluster medium. 

The timescale over which, following an instantaneous burst 
of SF, half of the Fe is released to the interstellar medium 
varies between ~ 0.3 to 3 Gyr for the wide variety of SNIa 
progenitor models considered here. Accordingly, the formation 
timescale of systems which exhibit an enhancement of a ele- 
ments with respect to Fe is rather uncertain, and depends on 
the SNIa model. Quantitatively, the actual constraint on the 
formation timescales also varies with the duration of the star 
formation episode in the system (Matteucci & Recchi 2001 ). 
Preliminary computations show that, in a star forming system, 
such timescale may range between 1 and several Gyr. This 
problem will be discussed in a forthcoming paper. 

The evolution of the gas flows in Ellipticals depends on the 
balance between the rate of mass return and the SNIa rate, past 
a burst of SF. The former scales with time as oc r 13 , while, at 
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delay times greater than ~ 1 Gyr, the fi d functions scale as oc t s 

with s 1 for the DD WIDE, 1.2 for the DD CLOSE, 

~ —1.6 for the SD Sub-Chandra. Therefore, it appears that the 
secular evolution of the SNIa rate past a burst of SF is critically 
close to the evolution of the rate of mass return, and that the 
fate of the gas in Ellipticals is very sensitive to the SNIa pro- 
genitor's model. It is also worth noticing that the shape of the 
/i a function is different from a simple power law, as is adopted 
in Ciotti et al. J1991I to model the gas flows in ellipticals. In 
particular, the presence of the wide maximum phase at inter- 
mediate epochs will impact on the dynamical evolution of the 
gas. 

In this paper, the emphasis has been put on the intercompar- 
ison of the various models for the SNIa progenitor. Actually, 
all different channels could contribute to the SNIa events, each 
with its own probability, as in the realizations of the population 
synthesis models. Some diversities of the observational prop- 
erties of SNIa have been found in the literature, which support 
this notion (e.g. Branch 2005 Benetti et al. 2005 ). In particular 
the different luminosity at maximum, and the different decline 
rate of the light curve, as measured by the the Amis parameter 
of Phillips ( 1993 1, of the events in early and late type galaxies 
(Delia Valle & Panagia Il992l van den Bergh & Pazder [19921 
Hamuy et al. 119961 Garnavich & Gallagher I2005> could be 
related to different typical progenitors. If both the Single and 
Double Degenerate channels are at work with similar total re- 
alization probabilities, in early type galaxies the DD explosions 
should prevail over SD events, since the distribution function of 
the delay times of the latter declines fast. In late type galaxies, 
instead, all channels should contribute to the current rate, with 
a larger proportion of SNIa from the SD channel, due to their 
high rate at early epochs. The formalism presented in this paper 
allows a straightforward exploration of the effect of a mixture 
of progenitors, e.g. by varying the relative Ai a realization prob- 
abilities. Eventually, it will be possible to constrain the mixture 
of progenitors by modeling the evolution of the SNIa rate in 
galaxies of different types, and considering as well all the other 
consequences on the large scales. 
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Appendix A: Equations used to derive the f^ D (r) 
function 

A .1. An approximate relation for the gravitational delay 

The gravitational delay time is given by: 
0.1 5 A 4 



T 2W — 



-Gyr. 



(A.l) 



(mm + m 2 R)m m m 2 R 
The mass dependent term can be written as 

/ = M DD ' m 2R - m\ R ■ M D D 

with M DD =mi R +m2R, i.e. the total mass of the DD system. 
Given that WD masses range between 0.6 and 1.2 M , M DD 
goes from 1 .2 to 2.4 M Q . f is a family of parabolas, both viewed 
as a function of m 2 R, and of Mud- These families are shown in 
Fig. iA.l\ . Each parabola in panel (a) has an absolute max- 
imum in m 2R = 0.5M DD , of / max = 0.25 • M DD . This locus 
is drawn in panel (a) as the dashed line, and in panel (b) as 
the uppermost thick curve. Not all combinations of (Mdd,'«2r) 




2.6 



Fig. A.l. Mass dependent term of r gw plotted versus »J2r for 
selected values of Mud (panel a), and plotted versus Mdd for 
selected values of m 2 R (panel b). The dashed curve in panel 
(a) connects the maxima of the parabolas. Along this line 
'«2r='«ir- The shaded area delimits the parameter space of 
SNIa producers. 



are acceptable, and the shaded area in the two panels show 
the allowed parameter space. In addition to the vertical lim- 
its 1.4 < M DD /M < 2.4, and 0.6 < m 2R /M e < 1.2, it is further 
required that ;7Jir(= Mdd - mod ranges between 0.6 M and 
1.2 M . The more restrictive criterion »J2r < wjr < 1.2M is 
unnecessary. Notice that, occasionally, the evolution in Close 
Binaries produces mi R < m 2R (even though, by definition, 
m\ > m 2 ), e.g. in some cases of conservative RLO when the 
secondary may become more massive than the original pri- 
mary, thereafter leaving a more massive WD remnant. 

The loci wir=O.6M and miR=1.2M are shown in 
Fig. IA.1I along each parabola, only the portion included be- 
tween the two loci is acceptable. Basically, at given Mud we 
exclude those values of m 2R which imply miR< 0.6 M , which 
is a He WD; similarly, we exclude those values of m 2 R which 
imply mjR> 1.2 M Q , taken here as the upper limit to the mass 
of a WD. 

As clearly visible in panel b) of Fig. IA.ll the term / is much 
more sensitive to M DD rather than to m2R, and there is an almost 
one to one correspondence between Muu and / over the whole 
parameter space of double CO WDs. In the computations pre- 
sented here I approximate the function / with its maximum 
value of 0.25 • M DD , which is the upper envelope of the family 
of parabolas in panel (b). This approximation leads to Eq. (I22> . 
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A.2. Shrinkage following the RLO 

Given the uncertainty of the results of the Common Envelope 
evolution I have considered two alternatives. In one case, both 
mass transfer phases are regulated through the standard CE Eq. 
i26\. which results into: 



.if m d m 
— = 0.5a ce — t- 



A; 



+ 0.5 a ce m 



(A.2) 



where m d indicates the mass of the donor, A indicates the sep- 
aration of the system, and the subscripts i and / refer to before 
and after the RLO, respectively. Further, m is the mass of the 
companion, and is the Roche Lobe radius in units of the ini- 
tial separation, which is adopted from Eggleton (1998), and is 
generally ~ 0.4 - 0.5. The blue and cyan dots in Fig.|4]result 
from the twofold application of this equation, first with 

mf=fenin\, 111=1112, and m d =OTn? as given by Eg. (I19> : 
second with 

m f-fen nil, m=mm, and ra d =/«2R as given by Eq. ( I19> . 

In the alternative evolutionary scheme, the first mass trans- 
fer is considered parametrized by the Envelope Ejection rela- 
tion in Nelemans et al. (2001 1, reported here as Eq. J27> . Since, 
in the same notation adopted above: 



mfm 



GAi 



1 K 



+ m 



rrAm 



GA f 



f A d 



+ m 



AM = mf - mf 



it follows 



A\ 



(m d \ 


2 


mf 


+ m 






m d 

i 


+ m , 



111; + 111. 



m d -m d ^ 
1-y-i 



'f 



(A.3) 



Eq. \A3\ is applicable if 



y- 



< l 



otherwise the system would loose more than 100% of its orig- 
inal angular momentum. With a little algebra, the last relation 
becomes: 



m d ^ 



1 



that is a lower limit to the mass ratio of the system at RLO. 
Typically, the quantity (mf/mf) (which is the core mass fraction 
of the donor at RLO) ranges between 0.15 to 0.25 for mf be- 
tween 8 and 2 M B . Then, q 0.2,0.4,0.6 for y = 1.5, 1.75,2. 
Close to this lower limit, the Af/Ai ratio from Eq. iA.3\ be- 
comes very small, and indeed smaller than that given by Eq. 

The red and magenta points in Fig.0]are generated with the 
following prescriptions: similar to Nelemans et al. (J2001 1, the 
first mass transfer results into a shrinkage given by Eq. JA.2t 
if the mass ratio is smaller than q, or the maximum Af/Aj from 
Eq. iA.2i and Eq. iA.'H if the mass ratio is larger than q. At the 
second mass transfer, Eq. jA.2t is applied. 



The product of the two Af jA\ ratios resulting from the first 
and the second mass transfer phases naturally equals the ratio 
between the final and the original separation of the close binary 
(A/A ) plotted in Fig.|H 

A.3. Function g(r, r n ) for WIDE DDs 

In this section I derive an expression for the fraction of systems 
which, having a nuclear delay r n , have a total delay smaller 
than r for the WIDE DD evolutionary scenario. Under the as- 
sumptions that the total binary mass (Mud) and the separation 
of the DD system (A) are independent variables, the contribu- 
tion to gravitational delay r gw from systems with separations in 
the range (A, A+dA) is: 



dn(-r gw ,A) = n(A) n(M DD ) dA dM DD 



(A.4) 



where n(A) and n(M DD ) are the distribution functions of the 
separations and total binary mass of the progeny of systems 
born with a secondary mass whose nuclear timescale is equal 



to T n , and the variables combine so that r RW =0.6- 



n(T gw )dT gw = dr g 



I n(A) n 
Ja„ 



(M dd ) 



8Mr 



gw 



dA 



, It follows: 



(A.5) 



where A n and A x define the range of separations which lead to 
the same gravitational delay r gw . Given the relation between 
gravitational delay, DD mass and separation, A n and A x are di- 
rectly related to the minimum and maximum M DD : 



A„ = 



M 0.75 0.25 
DD,n gw 

0.6 025 



M 0.75 0.25 
. _ DD,x gw 

x " 0.6 25 ' 



Fig. IA. 21 shows plausible limits for Mod, illustrating that 
the heavier n%2 is, the heavier its Mod progeny. To proceed, 
formulations for n(A) and «(Mdd) need to be specified: since 
the WIDE DD case is meant to describe a situation in which 
the evolution produces DDs in a wide range of separations, and 
since the relevant range of final separations is rather narrow 
(from 0.5 to 4.5 R e ), as a convenient parametrization I adopt: 



n(A) <x A A . 



(A.6) 



As for the distribution of the DD masses, I explore the conse- 
quences of two extreme assumptions: 

( 1 ) all values are equally probable: «(Mdd) = constant within 



M DD , n = max(1.4,m 2R + 0.6); M DDiX = m 2R + 1.2 
with ni2R given by Eq. dl9> : 

(2) the distribution is peaked around the value 10 (see Fig. lA.2> : 



M DD = 1.4 + (m 2 -2)/6. 

With respect to case (1), notice that, although one could enforce 
Mdd.ii = 2m2R, there are evolutionary paths which could lead 



10 This arbitrary relation is used just to explore the effect of a sys- 
tematic decrease of M DD as m 2 decreases. 
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2.5 



1 .5 




Fig. A.2. Limits for the total mass of the DD systems produced 
by primordial binaries with secondary mass mi, leaving a rem- 
nant m2R. The line labeled Mod is used in the computations 
to explore the effect of a tight correlation between mi and its 
progeny M DD . 

to ;«2r > miR, e.g. when after the 1st RLO the separation is 
larger than the primordial one. 

Case (2) corresponds to assuming that there is a tight corre- 
spondence between mi and the mass of the DD system, so that, 
at any mi (i.e. t„), Muu is specified, and the distribution of the 
gravitational delays of systems with given r n is: 

«(T gw )dT gw <x n(A)6A with A = (M dd )°- 75 t° w 25 /0.6 - 25 . 

With some algebra, Eq (IA.5> yields (for/? a i 1 — 1): 



"UgwJ ^ J 1,2 T gw 

where 

f W „ ^0.75+0.7% 
J\ DD 

for a narrow distribution of Muu around Mud, while 



J £ DD,x 



•Af, 



1 .75+0.75j8 a 
DD,n 



(A.7) 



(A.8) 



(A.9) 



for a wide distribution of Muu- 

Inserting Eq. &A.1\ in Eq. ( I25> . the g(r, r n ) function is de- 
rived as: 



g(r, r n ) 







(t-t.^-i* .) 

v u/ gw,i 



for T< T[ 

for Tl < T < T2 



(A. 10) 



T2 



t„ + r. 



gw.Xi 



with yS a = 0.25(1 +/3 a ), n = r n + r gw ,; 
and where the third branch has been neglected because it is 
assumed assume that the maximum gravitational delay is larger 
than the Hubble time for every r n . 
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Fig.A.3. Lower panel: factors in the RHS of Ea. JA.10l for a 
flat distribution of A. The dashed lines show the 7 dependent 
term and are labeled with the t value in Gyr. Notice how, as r 
increases toward t v =1 Gyr the baseline on r n increases, while 
for t > 1 the curves just shift upward. Upper panel: resulting g 
function in arbitrary units. Solid (dotted) lines result from using 
Eq.(E3DJwith/W(/W), 

The number of systems which have a nuclear delay r n and 
a total delay up to r is proportional to two factors: one (/*) de- 
scribing the systematics with Muu, the other scaling with the 
width of the parameter space in r gw . The dependence on the 
distribution of the separations appears in the exponents of both 
factors. The lower panel of Fig. IA.3l illustrates these two fac- 
tors for jS a =0, which corresponds to a flat distribution of final 
separations. The r dependent quantity 1 1 (dashed lines) is larger 
when the total delay is larger: more systems merge within a 
longer total delay. At fixed t, this factor decreases for increas- 
ing r n , as the available range in r gw decreases. At long total 



A small value of r gwi =0.001 Gyr has been used in this figure 
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delay times this effect becomes less important, since the upper 
limit to r n is of 1 Gyr only. The solid and the dotted lines in 
the lower panel of Fig. IA. 31 show respectively the and / 2 W 
factors, which account for the systematics of the gravitational 
delay with the nuclear delay: longer r n correspond to lower mi, 
and then to less massive Mod- At fixed A, systems with lower 
Mud are diluted over a larger r gw range: this effect is reflected 
on the decreasing f?t factors with increasing t„, and is more 
pronounced for because of the additional systematics with 
the range (Mdd.x - Muo,n), which gets smaller as r n increases 
(see Fig. 

The upper panel in Fig. IA.3I shows the resulting (non- 
normalized I2 ) function g(r, T n ) for the two considered distri- 
butions of Mod- By construction, the g functions are zero for 

Tn 



> t — r gW; i (see Eq. ( lA.lQt ). The fraction of systems within a 
given total delay time r decreases as r n increases, a dependence 
which is more pronounced when the variation in the range of 
Mud with is accounted for (i.e. when using / 2 W ). At any 
T n the fraction of systems with a delay time smaller than t in- 
creases with r. The variation of g with t„ and r depends on the 
distribution function of the final separations A, as discussed in 
Sect. 4. 3 and illustrated in Fig. [5] 



A.4. The differential distribution function of the delay 
times 

It its generic form, the Leibniz integral rule is: 

d C m f < w f b(z)d f 

— f(z,x)dx= — 



■ b(z) r b<z) df db da 

f(z, x)dx = \ -^-dx + f(z, b)— - f(z, a)—. 



a(z) dz dz 

This can be applied to Eq. i32i , with 

z = T X = T n f(z, x) = n(r n ) ■ g(j, T n ) 

a = r^j b = min(r x , t). 

Let's consider in turn the three additive terms at the right hand 
side of the Leibniz rule. The last term is equal to: 

dTni 

-["(Tn)g(T,T n )] Tn=Tn . — - =0 



because T n ,i is constant. 
The second term is equal to: 



[n(r a )g(T, T n )] Tn=min(T , ix>T; 



because 



dmin(r n , x ,T) 
dr 



(i) in t > r n x : = since t„ x is constant; 

(ii) in t < T n x : g(j, r n = t) = by construction, since the g 
function is null in r n > t - r gWi i (see Eq. d25l). 

Therefore, only the first term is left: 



-L 



min(r nx ,r) 



n(T n ) 7- dr n . 
or 



(A.ll) 



12 A proper definition requires ^(r gwx +t„, t„) = 1 at every r n , which 
is a relation between r gw x and r n . Since for the WIDE DD scheme 
we consider r < T gw-X for every r n the normalization factor is not 
important. 



The g(T, r n ) function is continuous (see Eq. \25i and thus its 
derivative can be computed in all its points. However, the func- 
tion presents cusps in r=r n +T gw ; and in T=T n +T gw x , where the 
derivative will not be continuous. 
Eq. iA~TTA ( WIDE DDs) can be written as: 







for T n > T - T„ W; j 



g(r, r n ) oc y ^ [fr _ _ for ^ ^ T _ ^ 

for delay times up to the Hubble time. Upon derivation this 
equation yields: 



dj_ K I % (T - rn)™ for r n j r - r gw , ; ^ 



forr n >T-T„ w ,i. 



Inserting Eq. jA.12l into Eq. flA.lll the distribution function 
of the delay times for the double degenerate systems in the 
WIDE DD evolutionary scheme is obtained as: 



Xmin(r nx ,r) 
«(T n )S W (T,T n )dT n 



(A.13) 



with 



S w (r,r n ) = 



/^(T-Tj- ™^ forTn < T _ TgWii 







for r n > r - T RW ,i 



For the CLOSE DD formulation, it's convenient to notice that 
the third branch in Eq. (13 1 i requires t > r n + r gw , x (T n ), and 
therefore exists only if the total delay time considered is suffi- 
ciently long: t > Tnj + T gw x (r n i). Therefore, I re-write Eq. i3l\ 

as: 

- if r < r n ,i + T gWjX (r„a): 



g(T, T n ) 







- if t > T n ,i + T n , x (T n ,i): 



g(j, T n ) 



1 







' gw,x ' gw ; 



for T n < T - T gw ,i 
for T n > T - T gw ,i 



for T„ < T* 



for T* < T n < T - T gw ,i 
for T n > T - T gw ,i 



(A. 14) 



(A.15) 



where t* is the solution of the equation t = t„+ T gw>x (T n ). The 
derivative is then: 

- if T < T^j + T gWjX (T n ,i): 



- if t > T n ,i + T„ iX (T n ,i): 



^ for T n < T - T gWji 
forTn > T - T gw i. 



dT 







for T n < T* 

(T-Tnf* f OI - T * < T < T _ T . 

AUA *n — 'n — ' ' gw,i 



for t„ > t - T gWil 

The last two equations can be written in a compact form as: 



dr 



/ C (T, T„) for T n < T - T gw>1 





for T n > T - T„ 
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with 



r(r,T n ) = 



for r < r n ,i + r gw>x (T n ,i): 



!+/%_ l*fit 
' ' ' * gw,i 



forT„<r* 

for T > Tnj + T gw>x (T n ,i): <! (t-tJh f * 



Inserting the last two equations into Eq. l lA.ll> . and split- 
ting the integration over r n in two parts, one from T n> i to t* and 
the other from r* to the upper limit (min(r n x ,T)), the distribu- 
tion of the delay times is obtained as: 



min(r nx ,r) 



n(j n )S (T,T n )dr n 



(A. 16) 



with 



5 L (T,T n )= r^-r, 



(t-tJ" for T < T _ T . 

+/)» l+^g lul 'n ' 'gw,i 







for T n > T - T gw ,i 



and 

^~n,inf 



T n>i for T < Tjy + T gw , x (T n ,i) 

t* for t > T nji + r gw>x (T n>i ) 



(A. 17) 



(A.18) 



